Parratt
struct, parse the XML file, and if needed, load the data.
Functions | |
Parratt * | initParratt (char *filename) |
Returns a newly allocated Parratt struct. | |
void | initData (struct _data *data) |
Loads the data from file. | |
int | initMethod (char *type, char *algo) |
Initializes the fit method types. | |
void | initResolution (struct _other *other) |
Loads the resolution values from file. |
Parratt* initParratt | ( | char * | filename | ) |
Initializes a newly allocated Parratt struct. The struct is loaded with the XML input by calling loadXML()
. If the do_fit flag is "yes", then the data to be fit is loaded as well by calling initData()
.
00027 { 00028 unsigned int i, j; 00029 Parratt *par; 00030 00031 /*--------------------------------- 00032 * Allocate a new Parratt struct. 00033 */ 00034 par = (Parratt *)malloc(sizeof(Parratt)); 00035 00036 /*--------------------------------- 00037 * Set the misc variables to defaults. 00038 */ 00039 par->s.doit = 0; 00040 par->s.qmin = 0.0; 00041 par->s.qmax = 0.3; 00042 par->s.numq = 200; 00043 par->s.numz = 200; 00044 for (j=0; j<64; j++) 00045 { 00046 par->s.refl[j] = 0; 00047 par->s.prof[j] = 0; 00048 } 00049 sprintf(par->s.refl, "reflectivity.txt"); 00050 sprintf(par->s.prof, "profile.txt"); 00051 par->s.pntp = PARRATT; 00052 00053 /*--------------------------------- 00054 * Zero the layers struct 00055 */ 00056 par->l.num = 0; 00057 for (i=0; i<MAX_LAYER; i++) 00058 { 00059 for (j=0; j<64; j++) 00060 { 00061 par->l.nam[i][j] = 0; 00062 } 00063 par->l.idx[i] = 0; 00064 par->l.typ[i] = 0; 00065 } 00066 00067 /*--------------------------------- 00068 * Zero the model struct 00069 */ 00070 par->m.num = 2; 00071 par->m.idx = 0; 00072 for (i=0; i<MAX_PARAM; i++) 00073 { 00074 par->m.par[i] = 0.0; 00075 par->m.fix[i] = 1; 00076 } 00077 par->m.par[1] = 1.0; 00078 00079 /*--------------------------------- 00080 * Zero the other struct 00081 */ 00082 par->g.plz = POLARIZED_UP; 00083 par->g.rtp = CONSTANT; 00084 par->g.rsn = 0.0; 00085 par->g.rl = SINGLE; 00086 par->g.rsNum = 0; 00087 for (i=0; i<64; i++) par->g.rnam[i] = 0; 00088 00089 /*--------------------------------- 00090 * Zero the data structs 00091 */ 00092 for (i=0; i<64; i++) par->d1.nam[i] = 0; 00093 par->d1.idx = 0; 00094 par->d1.num = 0; 00095 par->d1.wgt = WEIGHT_NONE; 00096 00097 for (i=0; i<64; i++) par->d2.nam[i] = 0; 00098 par->d2.idx = 0; 00099 par->d2.num = 0; 00100 par->d2.wgt = WEIGHT_NONE; 00101 00102 /*--------------------------------- 00103 * Parse the XML file 00104 */ 00105 loadXML(filename, par); 00106 if (par->l.num < 2) error("Not enough layers"); 00107 00108 /*---------------------------------*/ 00109 /* 00110 * The 'do' attribute flags whether to load data and par. 00111 */ 00112 if (par->s.doit) 00113 { 00114 initData(&par->d1); 00115 if (par->s.pntp == POLARIZED) 00116 { 00117 initData(&par->d2); 00118 } 00119 } 00120 00121 /*---------------------------------*/ 00122 /* 00123 * Get resolution values from file, if required. 00124 */ 00125 if (par->g.rl==ARRAY) initResolution(&par->g); 00126 00127 return par; 00128 }
void initData | ( | struct _data * | data | ) |
This function loads the data from file. It is expecting two or three column data. Data order is unimportant.
00131 { 00132 unsigned int i, j; 00133 char buff[120]; 00134 unsigned int ndat; 00135 double dtmp; 00136 FILE *fp = NULL; 00137 00138 /*-------------------------------------------- 00139 * Open the data file 00140 */ 00141 fp = fopen(data->nam, "r"); 00142 if (fp==NULL) error("Cannot open data file"); 00143 00144 /*-------------------------------------------- 00145 * Count the number of data points 00146 */ 00147 ndat= 0; 00148 for (i=0; i<120; i++) buff[i] = 0; 00149 while(fgets(buff, 119, fp) != NULL) 00150 { 00151 j = sscanf(buff, "%lf %lf %lf", &dtmp, &dtmp, &dtmp); 00152 if ((j == 2) || (j == 3)) ndat++; 00153 for (i=0; i<120; i++) buff[i] = 0; 00154 } 00155 rewind(fp); 00156 00157 /*-------------------------------------------- 00158 * Allocate memory for the data 00159 */ 00160 data->num = ndat; 00161 data->xpt = (double *)calloc(data->num,sizeof(double)); 00162 data->ypt = (double *)calloc(data->num,sizeof(double)); 00163 data->sig = (double *)calloc(data->num,sizeof(double)); 00164 00165 /*-------------------------------------------- 00166 * Read in the data 00167 */ 00168 ndat= 0; 00169 for (i=0; i<120; i++) buff[i] = 0; 00170 while(fgets(buff, 119, fp) != NULL) 00171 { 00172 j = sscanf(buff, "%lf %lf %lf", &data->xpt[ndat], \ 00173 &data->ypt[ndat], \ 00174 &data->sig[ndat]); 00175 /* 00176 * If we want sig from file, but don't get it, 00177 * what should we do? For now, set sig=1. 00178 */ 00179 if ((j == 2) && (data->wgt == WEIGHT_ERR)) data->sig[ndat] = 1.0; 00180 if ((j == 2) || (j == 3)) ndat++; 00181 for (i=0; i<120; i++) buff[i] = 0; 00182 } 00183 fclose(fp); 00184 00185 /*-------------------------------------------- 00186 * Set the sig variables 00187 */ 00188 for (i=0; i<data->num; i++) 00189 { 00190 if (data->wgt == WEIGHT_NONE) data->sig[i] = 1.0; 00191 else if (data->wgt == WEIGHT_R) data->sig[i] = data->ypt[i]; 00192 else if (data->wgt == WEIGHT_SQRT) data->sig[i] = sqrt(data->ypt[i]); 00193 } 00194 00195 return; 00196 }
int initMethod | ( | char * | type, | |
char * | algo | |||
) |
Default fitting type is LEAST_SQUARES.
00199 { 00200 int i = 0, j = 0; 00201 /* 00202 * In case no functions were in the XML file, 00203 * set some defaults. 00204 */ 00205 if ((type==NULL) || (algo==NULL)) 00206 { 00207 fit_method.solv_type = &lstsq_main; 00208 fit_method.lstsq_type = gsl_multifit_fdfsolver_lmder; 00209 return GSL_SUCCESS; 00210 } 00211 /* 00212 * The solver and method types are declared here so that 00213 * the user could choose 00214 * "SORT OF" CHECKED FOR VALID ENTRIES 00215 */ 00216 if (strncmp(type,"MULTIDIM_MINIMIZAION", 20)==0) 00217 { 00218 fit_method.solv_type = &mdmin_main; 00219 i = 1; 00220 } 00221 else if (strncmp(type,"LEAST_SQUARES", 13)==0) 00222 { 00223 fit_method.solv_type = &lstsq_main; 00224 i = 2; 00225 } 00226 else if (strncmp(type,"SIMULATED_ANNEAL", 16)==0) 00227 { 00228 fit_method.solv_type = &siman_main; 00229 i = 3; 00230 } 00231 /* 00232 * Should we error out in this case? 00233 */ 00234 else 00235 { 00236 fit_method.solv_type = &lstsq_main; 00237 i = 2; 00238 } 00239 00240 /*------------------------------------------------------*/ 00241 if (strncmp(algo,"CONJUGATE_FR", 12)==0) 00242 { 00243 fit_method.mdmin_type = gsl_multimin_fdfminimizer_conjugate_fr; 00244 j = 1; 00245 } 00246 else if (strncmp(algo,"CONJUGATE_PR", 12)==0) 00247 { 00248 fit_method.mdmin_type = gsl_multimin_fdfminimizer_conjugate_pr; 00249 j = 1; 00250 } 00251 else if (strncmp(algo,"VECTOR_BFGS", 11)==0) 00252 { 00253 fit_method.mdmin_type = gsl_multimin_fdfminimizer_vector_bfgs; 00254 j = 1; 00255 } 00256 else if (strncmp(algo,"STEEP_DESC", 10)==0) 00257 { 00258 fit_method.mdmin_type = gsl_multimin_fdfminimizer_steepest_descent; 00259 j = 1; 00260 } 00261 else if (strncmp(algo,"LMSDR", 5)==0) 00262 { 00263 fit_method.lstsq_type = gsl_multifit_fdfsolver_lmsder; 00264 j = 2; 00265 } 00266 else if (strncmp(algo,"LMDR", 4)==0) 00267 { 00268 fit_method.lstsq_type = gsl_multifit_fdfsolver_lmder; 00269 j = 2; 00270 } 00271 /* 00272 * Should we error out in this case? 00273 */ 00274 else 00275 { 00276 fit_method.lstsq_type = gsl_multifit_fdfsolver_lmder; 00277 j = 2; 00278 } 00279 00280 /*------------------------------------------------------*/ 00281 if (((i==1) && (j==2)) || ((i==2) && (j==1))) 00282 { 00283 GSL_ERROR("Error: wrong method/algorithm pair.", GSL_EINVAL); 00284 return GSL_EINVAL; 00285 } 00286 00287 return GSL_SUCCESS; 00288 }
void initResolution | ( | struct _other * | other | ) |
This function loads the resolution data from the specified file. The function expects the data in two column format: the first is the q's, the second column the associated resolutions. The data MUST be in ascending order of q for the loading to work properly. However, not all q must be specified. The listing can be thought of as defining 'q-ranges,' with each number pair defining the resolution starting at that q, and ending at the next specified q.
00291 { 00292 unsigned int i, j; 00293 char buff[120]; 00294 unsigned int ndat; 00295 double dtmp; 00296 FILE *fp = NULL; 00297 00298 /*-------------------------------------------- 00299 * Open the resolution file 00300 */ 00301 fp = fopen(other->rnam, "r"); 00302 if (fp==NULL) error("Cannot open resolution file"); 00303 00304 /*-------------------------------------------- 00305 * Count the number of data points 00306 */ 00307 ndat= 0; 00308 for (i=0; i<120; i++) buff[i] = 0; 00309 while(fgets(buff, 119, fp) != NULL) 00310 { 00311 j = sscanf(buff, "%lf %lf %lf", &dtmp, &dtmp, &dtmp); 00312 if ((j == 2) || (j == 3)) ndat++; 00313 for (i=0; i<120; i++) buff[i] = 0; 00314 } 00315 rewind(fp); 00316 00317 /*-------------------------------------------- 00318 * Allocate memory for the data 00319 */ 00320 other->rsNum = ndat; 00321 other->rsArrayq = (double *)calloc(other->rsNum,sizeof(double)); 00322 other->rsArrayval = (double *)calloc(other->rsNum,sizeof(double)); 00323 00324 /*-------------------------------------------- 00325 * Read in the data 00326 */ 00327 ndat= 0; 00328 for (i=0; i<120; i++) buff[i] = 0; 00329 while(fgets(buff, 119, fp) != NULL) 00330 { 00331 j = sscanf(buff, "%lf %lf", &other->rsArrayq[ndat], \ 00332 &other->rsArrayval[ndat] ); 00333 if ((j == 2) || (j == 3)) ndat++; 00334 for (i=0; i<120; i++) buff[i] = 0; 00335 } 00336 fclose(fp); 00337 00338 00339 }