All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
smooth_gsl.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2009-2014, 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?)
50 
51  \future Possibly create a new gsl_bspline class which replaces
52  the GSL bspline workspace
53 
54  \future Allow user to probe chi squared and the covariance?
55  */
56  class smooth_gsl {
57 
58  public:
59 
60  smooth_gsl() {
62  }
63 
64  /// Begin using x-values from vector \c ix
65  smooth_gsl(const gsl_vector *ix) {
67  x=ix;
68  init();
69  }
70 
71  ~smooth_gsl() {
72  free();
73  }
74 
75  /// Set the number of coefficients
76  void set_ncoeff(int incoeffs) {
77  ncoeffs=incoeffs;
78  }
79 
80  /// Set order
81  void set_order(int order) {
82  norder=order;
83  }
84 
85  /// Set parameters
86  void set_pars(int incoeffs, int order) {
87  ncoeffs=incoeffs;
88  norder=order;
89  }
90 
91  /** \brief Set the x-values
92 
93  \comment
94  This is just a reminder to me that this function can't
95  be virtual since it is called in a constructor.
96  \endcomment
97  */
98  void set_x(const gsl_vector *ix) {
99  x=ix;
100  init();
101  }
102 
103  /** \brief Smooth data in \c y with errors \c e returning result \c ys
104  */
105  int smooth_data(const gsl_vector *y, const gsl_vector *e, gsl_vector *ys);
106 
107  /** \brief Smooth data in \c y returning result \c ys
108  */
109  int smooth_data(const gsl_vector *y, gsl_vector *ys);
110 
111  protected:
112 
113  /// Number of free coefficients for spline
114  size_t ncoeffs;
115 
116  /// Order of spline to be used (4=cubic)
117  size_t norder;
118 
119  /// internally calculated, number of "segment" to split the data into
120  size_t nbreak;
121 
122  /// True of the x values have been set
123  bool x_set;
124 
125  /// Spline workspace
126  gsl_bspline_workspace *bw;
127 
128  /// Spline temporary vector
129  gsl_vector *B;
130 
131  /// Parameters of linear fit, y=X*c
132  gsl_vector *c;
133 
134  /// Linear fit workspace
135  gsl_multifit_linear_workspace *mw;
136 
137  /// Workspace for spline fitting
138  gsl_matrix *X;
139 
140  /// Values of the independent variable
141  const gsl_vector *x;
142 
143  /// Covariance matrix
144  gsl_matrix *cov;
145 
146  /// Construct un-weighted fit
147  int fit(const gsl_vector *y);
148 
149  /// Construct weighted fit
150  int fit_errors(const gsl_vector *y, const gsl_vector *e);
151 
152  /// calculate smoothed curve value for a certain xi
153  double calc_for_x(double xi);
154 
155  /** \brief Allocate memory and initialize splines
156 
157  \comment
158  This is just a reminder to me that this function can't
159  be virtual since it is called in a constructor.
160  \endcomment
161  */
162  int init();
163 
164  /// Set default values and zero pointers
165  void init_pointers_and_defs();
166 
167  /// Free memory
168  int free();
169 
170  };
171 
172 #ifndef DOXYGEN_NO_O2NS
173 }
174 #endif
175 
176 #endif
double calc_for_x(double xi)
calculate smoothed curve value for a certain xi
void init_pointers_and_defs()
Set default values and zero pointers.
const gsl_vector * x
Values of the independent variable.
Definition: smooth_gsl.h:141
void set_order(int order)
Set order.
Definition: smooth_gsl.h:81
void set_x(const gsl_vector *ix)
Set the x-values.
Definition: smooth_gsl.h:98
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:138
int free()
Free memory.
size_t nbreak
internally calculated, number of "segment" to split the data into
Definition: smooth_gsl.h:120
void set_pars(int incoeffs, int order)
Set parameters.
Definition: smooth_gsl.h:86
gsl_matrix * cov
Covariance matrix.
Definition: smooth_gsl.h:144
size_t ncoeffs
Number of free coefficients for spline.
Definition: smooth_gsl.h:114
bool x_set
True of the x values have been set.
Definition: smooth_gsl.h:123
int fit(const gsl_vector *y)
Construct un-weighted fit.
gsl_bspline_workspace * bw
Spline workspace.
Definition: smooth_gsl.h:126
size_t norder
Order of spline to be used (4=cubic)
Definition: smooth_gsl.h:117
gsl_vector * B
Spline temporary vector.
Definition: smooth_gsl.h:129
void set_ncoeff(int incoeffs)
Set the number of coefficients.
Definition: smooth_gsl.h:76
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:56
smooth_gsl(const gsl_vector *ix)
Begin using x-values from vector ix.
Definition: smooth_gsl.h:65
gsl_vector * c
Parameters of linear fit, y=X*c.
Definition: smooth_gsl.h:132
gsl_multifit_linear_workspace * mw
Linear fit workspace.
Definition: smooth_gsl.h:135
int init()
Allocate memory and initialize splines.

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).
Hosted at Get Object-oriented Scientific Computing
Lib at SourceForge.net. Fast, secure and Free Open Source software
downloads..