00001 #ifndef RECIPES_H 00002 #define RECIPES_H 00003 /******************************************************************************* 00004 * E.S.O. - VLT project 00005 * 00006 * "@(#) $Id: recipes.h,v 1.9 2005/07/04 15:33:15 amodigli Exp $" 00007 * 00008 * who when what 00009 * -------- -------- ---------------------------------------------- 00010 * schreib 05/06/00 created 00011 */ 00012 00013 /************************************************************************ 00014 * recipes.h 00015 * some numerical recipes 00016 *---------------------------------------------------------------------- 00017 */ 00018 00019 /* 00020 * header files 00021 */ 00022 00023 /* 00024 #include <stdio.h> 00025 #include "spiffi_types.h" 00026 */ 00027 #include <qfits.h> 00028 #include <cpl.h> 00029 #include <inttypes.h> 00030 #include <float.h> 00031 #include <math.h> 00032 00033 #include "eclipse.h" 00034 00035 00036 /*---------------------------------------------------------------------------- 00037 * defines 00038 *--------------------------------------------------------------------------*/ 00039 /* definitions of initial values for lsqfit_c in linefit() (wave_calibration) */ 00040 #define XDIM 1 /* dimension of the x values */ 00041 #define TOL 0.001 /* fitting tolerance */ 00042 #define LAB 0.1 /* labda parameter */ 00043 #define ITS 200 /* maximum number of iterations */ 00044 #define MAXPAR 4 /* number of free parameters */ 00045 #define LABFAC 10.0 /* labda step factor */ 00046 #define LABMAX 1.0e+10 /* maximum value for labda */ 00047 #define LABMIN 1.0e-10 /* minimum value for labda */ 00048 00049 /*---------------------------------------------------------------------------- 00050 * Function ANSI C prototypes 00051 *--------------------------------------------------------------------------*/ 00052 00053 00054 00055 /*-------------------------------------------------------------------------*/ 00068 /*--------------------------------------------------------------------------*/ 00069 00070 OneCube * 00071 copy_cube_range( 00072 OneCube *src_cube, 00073 const int z_min, 00074 const int z_max 00075 ); 00076 00077 /*---------------------------------------------------------------------------- 00078 Function : myfit() 00079 In : set of data points, individual standard deviations, 00080 parameters of a straight line, chi-square, goodness-of- 00081 fit probability 00082 Out : void 00083 Job : fits a straight line to a set of x, y data points by 00084 minimizing chi-square. Numeric. Rec. page 665 00085 ---------------------------------------------------------------------------*/ 00086 00087 void myfit (float x[], float y[], int ndata, float sig[], int mwt, float *a, 00088 float *b, float *siga, float *sigb, float *chi2, float *q) ; 00089 00090 /*---------------------------------------------------------------------------- 00091 Function : median() 00092 In : array of pixelvalues starting at index 0, 00093 number of array elements 00094 Out : median of pixelvalues 00095 Job : sorts an array into ascending numerical order by 00096 using pixel_qsort(), and derives the median depending 00097 on an odd or even number of array elements 00098 ---------------------------------------------------------------------------*/ 00099 00100 pixelvalue median(pixelvalue * array, int n) ; 00101 00102 /*<python>*/ 00103 float f_median(float * array, int n) ; 00104 /*</python>*/ 00105 00106 /*---------------------------------------------------------------------------- 00107 Function : gaussian() 00108 In : position array xdat, parameter list parlist, number of 00109 parameters in the list npar 00110 The parameters are: 00111 parlist(0): Amplitude 00112 parlist(1): FWHM 00113 parlist(2): x0, center of gauss versus center of column 00114 parlist(3): Zero level 00115 Out : function value of a 1-d gaussian: 00116 F(x) = par(0) * EXP ( -4.0 * LN(2.0) * 00117 [((x - parlist(2))/parlist(1))^2]) + parlist(3) 00118 FWHM^2 = 8 * ln (2) * sigma^2 00119 Job : calculates the value of a gaussian with parameters 00120 parlist at the position xdat 00121 ---------------------------------------------------------------------------*/ 00122 00123 float gaussian ( float * xdat, float * parlist/*, int * npar*/ ) ; 00124 00125 /*---------------------------------------------------------------------------- 00126 Function : gaussian_deriv() 00127 In : position array xdat, parameter list parlist, number of 00128 parameters in the list npar 00129 The parameters are: 00130 parlist(0): Amplitude 00131 parlist(1): Axis(FWHM) 00132 parlist(2): x0, center of gauss versus center of column 00133 parlist(3): Zero level 00134 derivative value of gaussian at position xdat: dervs 00135 dervs[0]: partial derivative by the amplitude 00136 dervs[1]: partial derivative by FWHM 00137 dervs[2]: partial derivative by the x position (parlist[2]) 00138 dervs[3]: partial derivative by the zero level 00139 Out : nothing, void 00140 Job : calculates the partial derivatives for a gaussian with 00141 parameters parlist at position xdat 00142 ---------------------------------------------------------------------------*/ 00143 00144 void gaussian_deriv( float * xdat, float * parlist, float * dervs/*, int * npar*/ ) ; 00145 00146 /*---------------------------------------------------------------------------- 00147 Function : lsqfit_c() 00148 In : xdat: position, coordinates of data points. 00149 xdat is 2 dimensional: XDAT ( XDIM, NDAT ) 00150 xdim: dimension of fit 00151 ydat: data points 00152 wdat: weights for data points 00153 ndat: number of data points 00154 fpar: on input contains initial estimates of the 00155 parameters for non-linear fits, on output the 00156 fitted parameters. 00157 epar: contains estimates of the errors in fitted parameters 00158 mpar: logical mask telling which parameters are free (non-zero) 00159 and which parameters are fixed (0) 00160 npar: number of function parameters ( free + fixed ) 00161 tol: relative tolerance. lsqfit stops when successive iterations 00162 fail to produce a decrement in reduced chi-squared less 00163 than tol. If tol is less than the minimum tolerance 00164 possible, tol will be set to this value. This means 00165 that maximum accuracy can be obtained by setting 00166 tol = 0.0. 00167 its: maximum number of iterations 00168 lab: mixing parameter, lab determines the initial weight 00169 of steepest descent method relative to the Taylor method 00170 lab should be a small value (i.e. 0.01). lab can only 00171 be zero when the partial derivatives are independent 00172 of the parameters. In fact in this case lab should be 00173 exactly equal to zero. 00174 Out : returns number of iterations needed to achieve convergence 00175 according to tol. When this number is negative, the fitting 00176 was not continued because a fatal error occurred: 00177 -1 too many free parameters, maximum is 32 00178 -2 no free parameters 00179 -3 not enough degrees of freedom 00180 -4 maximum number of iterations too small to obtain 00181 a solution which satisfies tol. 00182 -5 diagonal of matrix contains elements which are zero 00183 -6 determinant of the coefficient matrix is zero 00184 -7 square root of a negative number 00185 Job : this is a routine for making a least-squares fit of a 00186 function to a set of data points. The method used is 00187 described in: Marquardt, J.Soc.Ind.Appl.Math. 11. 431 (1963). 00188 This method is a mixture of the steepest descent method 00189 and the Taylor method. 00190 ---------------------------------------------------------------------------*/ 00191 00192 int lsqfit_c ( float * xdat, 00193 int * xdim, 00194 float * ydat, 00195 float * wdat, 00196 int * ndat, 00197 float * fpar, 00198 float * epar, 00199 int * mpar, 00200 int * npar, 00201 float * tol , 00202 int * its , 00203 float * lab ) ; 00204 00205 /*---------------------------------------------------------------------------- 00206 Function : nint() 00207 In : x : specified real value 00208 Out : nearest integer to specified real value 00209 Job : this routine determines the nearest integer to a 00210 specified real value 00211 ---------------------------------------------------------------------------*/ 00212 00213 int nint ( double x ) ; 00214 00215 /*---------------------------------------------------------------------------- 00216 Function : correlation() 00217 In : data1: first data array 00218 data2: second data array 00219 ndata: number of data elements in the arrays, 00220 both arrays must have the same number of 00221 elements. 00222 Out : integer shift of the second array with respect to 00223 the first array. 00224 2^32/2 if something went wrong. 00225 Job : this routine computes the cross correlation between 00226 two data arrays of the same size and returns the 00227 integer shift between both arrays 00228 ---------------------------------------------------------------------------*/ 00229 00230 int correlation ( float * data1, float * data2, int ndata ) ; 00231 00232 /*---------------------------------------------------------------------------- 00233 Function : clean_mean() 00234 In : array: data array to average 00235 n_elements: number of elements of the data array 00236 throwaway_low: percentage of low value elements to be thrown away before 00237 averaging 00238 throwaway_high: percentage of high value elements to be thrown away before 00239 averaging 00240 Out : the clean mean of a data array 00241 FLT_MAX in case of error 00242 Job : this routine computes the clean mean of a given data 00243 array that means the array is first sorted and 00244 a given percentage of the lowest and the highest values 00245 is not considered for averaging 00246 ---------------------------------------------------------------------------*/ 00247 /*<python>*/ 00248 float clean_mean( float * array, 00249 int n_elements, 00250 float throwaway_low, 00251 float throwaway_high ) ; 00252 /*</python>*/ 00253 00254 /*---------------------------------------------------------------------------- 00255 Function : errorPrint() 00256 In : text: text for e_error to print out 00257 Out : nothing 00258 Job : this routine gives python the chance to print out 00259 error messages 00260 ---------------------------------------------------------------------------*/ 00261 /*<python>*/ 00262 void errorPrint(char * text) ; 00263 /*</python>*/ 00264 00265 /*---------------------------------------------------------------------------- 00266 Function : commentPrint() 00267 In : text: text for e_comment to print out 00268 level: comment indentation level 00269 Out : nothing 00270 Job : this routine gives python the chance to print out 00271 comment messages 00272 ---------------------------------------------------------------------------*/ 00273 /*<python>*/ 00274 void commentPrint(char * text) ; 00275 /*</python>*/ 00276 00277 /*---------------------------------------------------------------------------- 00278 Function : convertZEROsTo0ForImages() 00279 In : input image 00280 Out : nothing 00281 Job : this function converts any ZEROs in an image to 0 00282 ---------------------------------------------------------------------------*/ 00283 /*<python>*/ 00284 void convertZEROsTo0ForImages(OneImage * im) ; 00285 /*</python>*/ 00286 00287 00288 00289 /*---------------------------------------------------------------------------- 00290 Function : convertZEROsTo0ForCubesRange() 00291 In : input cube 00292 Out : nothing 00293 Job : this function converts any ZEROs in a cube to 0 00294 ---------------------------------------------------------------------------*/ 00295 void convertZEROsTo0ForCubesRange(OneCube * cube,const int z_min,const int z_max); 00296 00297 00298 /*---------------------------------------------------------------------------- 00299 Function : convert0ToZEROsForCubes() 00300 In : input cube 00301 z_min min z to be considered 00302 z_max max z to be considered 00303 Out : nothing 00304 Job : this function converts any 0 in a cube to ZERO 00305 ---------------------------------------------------------------------------*/ 00306 void convert0ToZEROForCubesRange(OneCube * cube,const int z_min,const int z_max); 00307 /*---------------------------------------------------------------------------- 00308 Function : convertZEROsTo0ForCubes() 00309 In : input cube 00310 Out : nothing 00311 Job : this function converts any ZEROs in a cube to 0 00312 ---------------------------------------------------------------------------*/ 00313 /*<python>*/ 00314 void convertZEROsTo0ForCubes(OneCube * cube) ; 00315 /*</python>*/ 00316 00317 /*---------------------------------------------------------------------------- 00318 Function : convert0ToZEROsForImages() 00319 In : input image 00320 Out : nothing 00321 Job : this function converts any 0 in an image to ZEROs 00322 ---------------------------------------------------------------------------*/ 00323 void convert0ToZEROsForImages(OneImage * im) ; 00324 00325 /*---------------------------------------------------------------------------- 00326 Function : convert0ToZEROsForCubes() 00327 In : input cube 00328 Out : nothing 00329 Job : this function converts any 0 in a cube to ZERO 00330 ---------------------------------------------------------------------------*/ 00331 void convert0ToZEROForCubes(OneCube * cube) ; 00332 00333 /*---------------------------------------------------------------------------- 00334 Function : invertImage() 00335 In : input image 00336 Out : inverted input image (overwritten) 00337 Job : this function multiplies each image value by -1. 00338 ---------------------------------------------------------------------------*/ 00339 /*<python>*/ 00340 void invert(OneImage * im) ; 00341 /*</python>*/ 00342 00343 /*---------------------------------------------------------------------------- 00344 Function : xcorrel () 00345 In : line_i: the reference signal 00346 width_i: number of samples in reference signal 00347 line_t: candidate signal to compare 00348 width_t: number of samples in candidate signal 00349 half_search: half-size of the search domain 00350 delta: output correlation offset. 00351 maxpos: position index of the maximum of the output array 00352 xcorr_max: maximum value of the output array 00353 Out : correlation function, must be freed. 00354 Job : simple cross-correlation of two 1d-signals without FFT 00355 ---------------------------------------------------------------------------*/ 00356 00357 double * xcorrel( 00358 float * line_i, 00359 int width_i, 00360 float * line_t, 00361 int width_t, 00362 int half_search, 00363 int * delta, 00364 int * maxpos, 00365 double * xcorr_max 00366 00367 ) ; 00368 00369 float nev_ille(float [], float [], int, float, int *); 00370 00371 #endif 00373 /*--------------------------------------------------------------------------*/ 00374
1.2.13.1 written by Dimitri van Heesch,
© 1997-2001