smooth_gsl.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2009-2017, Marco Cammarata and Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 #ifndef SMOOTH_GSL_H
24 #define SMOOTH_GSL_H
25 
26 /** \file smooth_gsl.h
27  \brief File defining \ref o2scl::smooth_gsl
28 */
29 
30 #include <iostream>
31 #include <fstream>
32 #include <string>
33 #include <cmath>
34 #include <sstream>
35 
36 #include <gsl/gsl_bspline.h>
37 #include <gsl/gsl_multifit.h>
38 
39 #ifndef DOXYGEN_NO_O2NS
40 namespace o2scl {
41 #endif
42 
43  /** \brief Smooth a GSL vector using GSL bsplines
44 
45  \todo Needs a bit more error checking and
46  more documentation.
47 
48  \future Generalize to generic vector types. (Does this require
49  reworking the GSL linear fitting routines? Doesn't matter now,
50  the GSL linear fitting routines are now reworked.)
51 
52  \future Possibly create a new gsl_bspline class which replaces
53  the GSL bspline workspace
54 
55  \future Allow user to probe chi squared and the covariance?
56  */
57  class smooth_gsl {
58 
59  public:
60 
61  smooth_gsl() {
63  }
64 
65  /// Begin using x-values from vector \c ix
66  smooth_gsl(const gsl_vector *ix) {
68  x=ix;
69  init();
70  }
71 
72  ~smooth_gsl() {
73  free();
74  }
75 
76  /// Set the number of coefficients
77  void set_ncoeff(int incoeffs) {
78  ncoeffs=incoeffs;
79  }
80 
81  /// Set order
82  void set_order(int order) {
83  norder=order;
84  }
85 
86  /// Set parameters
87  void set_pars(int incoeffs, int order) {
88  ncoeffs=incoeffs;
89  norder=order;
90  }
91 
92  /** \brief Set the x-values
93 
94  \comment
95  This is just a reminder to me that this function can't
96  be virtual since it is called in a constructor.
97  \endcomment
98  */
99  void set_x(const gsl_vector *ix) {
100  x=ix;
101  init();
102  }
103 
104  /** \brief Smooth data in \c y with errors \c e returning result \c ys
105  */
106  int smooth_data(const gsl_vector *y, const gsl_vector *e, gsl_vector *ys);
107 
108  /** \brief Smooth data in \c y returning result \c ys
109  */
110  int smooth_data(const gsl_vector *y, gsl_vector *ys);
111 
112  protected:
113 
114  /// Number of free coefficients for spline
115  size_t ncoeffs;
116 
117  /// Order of spline to be used (4=cubic)
118  size_t norder;
119 
120  /// internally calculated, number of "segment" to split the data into
121  size_t nbreak;
122 
123  /// True of the x values have been set
124  bool x_set;
125 
126  /// Spline workspace
127  gsl_bspline_workspace *bw;
128 
129  /// Spline temporary vector
130  gsl_vector *B;
131 
132  /// Parameters of linear fit, y=X*c
133  gsl_vector *c;
134 
135  /// Linear fit workspace
136  gsl_multifit_linear_workspace *mw;
137 
138  /// Workspace for spline fitting
139  gsl_matrix *X;
140 
141  /// Values of the independent variable
142  const gsl_vector *x;
143 
144  /// Covariance matrix
145  gsl_matrix *cov;
146 
147  /// Construct un-weighted fit
148  int fit(const gsl_vector *y);
149 
150  /// Construct weighted fit
151  int fit_errors(const gsl_vector *y, const gsl_vector *e);
152 
153  /// calculate smoothed curve value for a certain xi
154  double calc_for_x(double xi);
155 
156  /** \brief Allocate memory and initialize splines
157 
158  \comment
159  This is just a reminder to me that this function can't
160  be virtual since it is called in a constructor.
161  \endcomment
162  */
163  int init();
164 
165  /// Set default values and zero pointers
166  void init_pointers_and_defs();
167 
168  /// Free memory
169  int free();
170 
171  };
172 
173 #ifndef DOXYGEN_NO_O2NS
174 }
175 #endif
176 
177 #endif
double calc_for_x(double xi)
calculate smoothed curve value for a certain xi
The main O<span style=&#39;position: relative; top: 0.3em; font-size: 0.8em&#39;>2</span>scl O$_2$scl names...
Definition: anneal.h:42
void init_pointers_and_defs()
Set default values and zero pointers.
const gsl_vector * x
Values of the independent variable.
Definition: smooth_gsl.h:142
void set_order(int order)
Set order.
Definition: smooth_gsl.h:82
void set_x(const gsl_vector *ix)
Set the x-values.
Definition: smooth_gsl.h:99
int smooth_data(const gsl_vector *y, const gsl_vector *e, gsl_vector *ys)
Smooth data in y with errors e returning result ys.
gsl_matrix * X
Workspace for spline fitting.
Definition: smooth_gsl.h:139
int free()
Free memory.
size_t nbreak
internally calculated, number of "segment" to split the data into
Definition: smooth_gsl.h:121
void set_pars(int incoeffs, int order)
Set parameters.
Definition: smooth_gsl.h:87
gsl_matrix * cov
Covariance matrix.
Definition: smooth_gsl.h:145
size_t ncoeffs
Number of free coefficients for spline.
Definition: smooth_gsl.h:115
bool x_set
True of the x values have been set.
Definition: smooth_gsl.h:124
int fit(const gsl_vector *y)
Construct un-weighted fit.
gsl_bspline_workspace * bw
Spline workspace.
Definition: smooth_gsl.h:127
size_t norder
Order of spline to be used (4=cubic)
Definition: smooth_gsl.h:118
gsl_vector * B
Spline temporary vector.
Definition: smooth_gsl.h:130
void set_ncoeff(int incoeffs)
Set the number of coefficients.
Definition: smooth_gsl.h:77
int fit_errors(const gsl_vector *y, const gsl_vector *e)
Construct weighted fit.
Smooth a GSL vector using GSL bsplines.
Definition: smooth_gsl.h:57
smooth_gsl(const gsl_vector *ix)
Begin using x-values from vector ix.
Definition: smooth_gsl.h:66
gsl_vector * c
Parameters of linear fit, y=X*c.
Definition: smooth_gsl.h:133
gsl_multifit_linear_workspace * mw
Linear fit workspace.
Definition: smooth_gsl.h:136
int init()
Allocate memory and initialize splines.

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).