00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 #ifdef HAVE_CONFIG_H
00027 #include <config.h>
00028 #endif
00029
00030
00039 #include <xsh_drl.h>
00040 #include <xsh_pfits.h>
00041 #include <xsh_pfits_qc.h>
00042 #include <xsh_utils.h>
00043 #include <xsh_data_order.h>
00044 #include <xsh_error.h>
00045 #include <xsh_utils.h>
00046 #include <xsh_msg.h>
00047 #include <xsh_data_pre.h>
00048 #include <xsh_data_instrument.h>
00049 #include <xsh_data_order.h>
00050 #include <xsh_data_wavesol.h>
00051 #include <xsh_data_resid_tab.h>
00052 #include <xsh_data_wavemap.h>
00053 #include <xsh_data_spectralformat.h>
00054 #include <xsh_model_io.h>
00055 #include <xsh_model_kernel.h>
00056 #include <cpl.h>
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078 void xsh_create_map( cpl_frame *dispsol_frame, cpl_frame *ordertab_frame,
00079 cpl_frame *pre_frame, xsh_instrument *instrument, cpl_frame **wavemap_frame,
00080 cpl_frame **slitmap_frame,const char* rec_prefix)
00081 {
00082 xsh_dispersol_list *dispersol_tab = NULL;
00083 xsh_pre *pre = NULL;
00084 char wavemap_tag[256];
00085 char slitmap_tag[256];
00086
00087 XSH_ASSURE_NOT_NULL( dispsol_frame);
00088 XSH_ASSURE_NOT_NULL( ordertab_frame);
00089 XSH_ASSURE_NOT_NULL( pre_frame);
00090 XSH_ASSURE_NOT_NULL( instrument);
00091 XSH_ASSURE_NOT_NULL( wavemap_frame);
00092 XSH_ASSURE_NOT_NULL( slitmap_frame);
00093
00094
00095 check( pre = xsh_pre_load( pre_frame, instrument));
00096 check( dispersol_tab = xsh_dispersol_list_load( dispsol_frame,
00097 instrument));
00098 sprintf(wavemap_tag,"%s_%s",
00099 rec_prefix,XSH_GET_TAG_FROM_ARM( XSH_WAVE_MAP_POLY, instrument));
00100 sprintf(slitmap_tag,"%s_%s",
00101 rec_prefix,XSH_GET_TAG_FROM_ARM( XSH_SLIT_MAP_POLY, instrument));
00102
00103 check( *wavemap_frame = xsh_dispersol_list_to_wavemap( dispersol_tab,
00104 ordertab_frame, pre, instrument,wavemap_tag));
00105 check( *slitmap_frame = xsh_dispersol_list_to_slitmap( dispersol_tab,
00106 ordertab_frame,
00107 pre, instrument,
00108 slitmap_tag));
00109
00110 cleanup:
00111 xsh_dispersol_list_free( &dispersol_tab);
00112 xsh_pre_free( &pre);
00113 return;
00114 }
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127 void xsh_create_model_map( cpl_frame* model_frame, xsh_instrument* instrument,
00128 const char* wtag, const char* stag,
00129 cpl_frame **wavemap_frame,
00130 cpl_frame **slitmap_frame,const int save_tmp)
00131 {
00132 xsh_xs_3 model_config;
00133
00134
00135 XSH_ASSURE_NOT_NULL_MSG( model_frame,"If model-scenario is 0 make sure that the input model cfg has at least one parameter with Compute_Flag set to 1 and High_Limit>Low_limit");
00136 XSH_ASSURE_NOT_NULL( instrument);
00137 XSH_ASSURE_NOT_NULL( wavemap_frame);
00138 XSH_ASSURE_NOT_NULL( slitmap_frame);
00139 XSH_ASSURE_NOT_NULL( wtag);
00140 XSH_ASSURE_NOT_NULL( stag);
00141
00142 check( xsh_model_config_load_best( model_frame, &model_config));
00143
00144
00145
00146
00147
00148
00149 check( xsh_model_binxy( &model_config, instrument->binx,
00150 instrument->biny));
00151
00152 check(xsh_model_maps_create(&model_config,instrument,wtag,stag,
00153 wavemap_frame,slitmap_frame,save_tmp));
00154
00155
00156
00157
00158
00159
00160
00161 cleanup:
00162
00163
00164 return;
00165 }
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194
00195 cpl_frame *
00196 xsh_create_poly_wavemap( cpl_frame *pre_frame,
00197 cpl_frame *wave_tab_2d_frame,
00198 cpl_frame *order_tab_frame,
00199 cpl_frame *spectral_format_frame,
00200 xsh_dispersol_param *dispsol_par,
00201 xsh_instrument * instrument,
00202 const char * wm_tag,
00203 cpl_frame **dispersol_frame,
00204 cpl_frame** slitmap_frame)
00205 {
00206 xsh_pre *pre = NULL;
00207 xsh_spectralformat_list *spec_list = NULL;
00208 xsh_wavesol *wave_tab_2d = NULL;
00209 cpl_frame *result = NULL ;
00210 float slit_step = 1.5;
00211 float lambda_step = 0.1;
00212
00213
00214 float sol_min_slit, sol_max_slit;
00215
00216 int i, idx, size, slit_size;
00217 double j, k;
00218 double *vlambda = NULL, *vslit = NULL, *vorder = NULL;
00219 double *pos_x = NULL, *pos_y = NULL;
00220
00221 xsh_wavemap_list *wm_list = NULL ;
00222 int binx;
00223
00224 char wm_name[256];
00225
00226
00227 XSH_ASSURE_NOT_NULL( wave_tab_2d_frame);
00228 XSH_ASSURE_NOT_NULL( order_tab_frame);
00229 XSH_ASSURE_NOT_NULL( spectral_format_frame);
00230 XSH_ASSURE_NOT_NULL( dispsol_par);
00231 XSH_ASSURE_NOT_NULL( instrument);
00232 XSH_ASSURE_NOT_NULL( dispersol_frame);
00233
00234
00235 if (pre_frame != NULL){
00236 check( pre = xsh_pre_load( pre_frame, instrument));
00237 }
00238
00239 check( binx = xsh_instrument_get_binx( instrument )) ;
00240
00241
00242 check( spec_list = xsh_spectralformat_list_load( spectral_format_frame,
00243 instrument));
00244 check( wave_tab_2d = xsh_wavesol_load( wave_tab_2d_frame,
00245 instrument ));
00246
00247 check( xsh_wavesol_set_bin_x( wave_tab_2d, binx ) ) ;
00248 check( xsh_wavesol_set_bin_y( wave_tab_2d, binx ) ) ;
00249
00250
00251
00252
00253
00254 sol_min_slit = wave_tab_2d->min_slit;
00255 sol_max_slit = wave_tab_2d->max_slit;
00256
00257
00258
00259
00260
00261
00262 slit_size = (int)((sol_max_slit-sol_min_slit)/slit_step)+1;
00263
00264 size=0;
00265 for( i=0; i< spec_list->size; i++){
00266 float lambda_min, lambda_max;
00267 int lambda_size;
00268
00269 lambda_min = spec_list->list[i].lambda_min_full;
00270 lambda_max = spec_list->list[i].lambda_max_full;
00271 XSH_ASSURE_NOT_ILLEGAL_MSG(lambda_max>= lambda_min,
00272 "lambda_max< lambda_min!! Check input spectralformat table, columns WLMAXFUL and WLMINFUL");
00273 lambda_size = (int)((lambda_max-lambda_min)/lambda_step)+1;
00274 size += lambda_size*slit_size;
00275 }
00276 xsh_msg_dbg_medium( "size %d", size ) ;
00277
00278 idx = 0;
00279 XSH_MALLOC( vorder, double, size);
00280 XSH_MALLOC( vlambda, double, size);
00281 XSH_MALLOC( vslit, double, size);
00282 XSH_MALLOC( pos_x, double, size);
00283 XSH_MALLOC( pos_y, double, size);
00284
00285 for( i=0; i< spec_list->size; i++){
00286 double absorder;
00287 float lambda_min, lambda_max;
00288
00289 absorder= (double)spec_list->list[i].absorder;
00290 lambda_min = spec_list->list[i].lambda_min_full;
00291 lambda_max = spec_list->list[i].lambda_max_full;
00292 xsh_msg("order %f lambda %f-%f",absorder, lambda_min, lambda_max);
00293
00294 for(j=lambda_min; j <= lambda_max; j+=lambda_step){
00295 for( k=sol_min_slit; k <= sol_max_slit; k+=slit_step){
00296 double x, y;
00297
00298 check( x = xsh_wavesol_eval_polx( wave_tab_2d, j, absorder, k));
00299 check( y = xsh_wavesol_eval_poly( wave_tab_2d, j, absorder, k));
00300 vorder[idx] = absorder;
00301 vlambda[idx] = j;
00302 vslit[idx] = k;
00303 pos_x[idx] = x;
00304 pos_y[idx] = y;
00305 idx++;
00306 }
00307 }
00308 xsh_msg_dbg_medium( "i %d idx %d / %d", i, idx, size);
00309 }
00310
00311 if ( xsh_debug_level_get() >= XSH_DEBUG_LEVEL_MEDIUM){
00312 FILE *test = NULL;
00313 int itest;
00314
00315 test = fopen( "wavemap_grid.log", "w");
00316 for( itest=0; itest < size; itest++){
00317 fprintf(test, "%f %f\n", pos_x[itest], pos_y[itest]);
00318 }
00319 fclose(test);
00320 }
00321
00322
00323 check( wm_list = xsh_wavemap_list_create( instrument));
00324 #if 0
00325 check( xsh_wavemap_list_compute( vlambda, pos_x, pos_y, size,
00326 vorder, wm_par, wm_list));
00327 #else
00328
00329
00330
00331
00332
00333
00334 check( xsh_wavemap_list_compute_poly( vlambda, vslit, pos_x, pos_y, idx,
00335 vorder, dispsol_par, wm_list));
00336 #endif
00337 sprintf(wm_name,"%s.fits",wm_tag);
00338 check( result = xsh_wavemap_list_save_poly( wm_list, order_tab_frame,
00339 pre, instrument, wm_tag, dispersol_frame, slitmap_frame));
00340 if ( pre != NULL){
00341
00342 check( cpl_frame_set_tag( result,wm_tag));
00343 check( cpl_frame_set_tag( *slitmap_frame,
00344 XSH_GET_TAG_FROM_ARM( XSH_SLIT_MAP_POLY, instrument)));
00345 }
00346 cleanup:
00347
00348 XSH_FREE( vorder);
00349 XSH_FREE( vlambda);
00350 XSH_FREE( vslit);
00351 XSH_FREE( pos_x);
00352 XSH_FREE( pos_y);
00353 xsh_pre_free( &pre);
00354 xsh_spectralformat_list_free( &spec_list);
00355 xsh_wavesol_free( &wave_tab_2d);
00356 xsh_wavemap_list_free( &wm_list);
00357 return result;
00358 }
00359
00360
00361
00362 cpl_frame*
00363 xsh_create_dispersol_physmod(cpl_frame *pre_frame,
00364 cpl_frame *order_tab_frame,
00365 cpl_frame* mod_cfg_frame,
00366 cpl_frame* wave_map_frame,
00367 cpl_frame* slit_map_frame,
00368 xsh_dispersol_param *dispsol_param,
00369 cpl_frame* spectral_format_frame,
00370 xsh_instrument* instrument,
00371 const int clean_tmp)
00372 {
00373
00374
00375 xsh_pre *pre = NULL;
00376 cpl_frame* disp_frame=NULL;
00377 xsh_spectralformat_list *spec_list = NULL;
00378
00379
00380 xsh_xs_3 model_config;
00381 int slit_size=0;
00382 int slit_step=1.5;
00383
00384
00385
00386 float sol_min_slit=0, sol_max_slit=0;
00387 int sol_min_order=0, sol_max_order=0;
00388 cpl_image* wmap_ima=NULL;
00389 cpl_image* smap_ima=NULL;
00390 cpl_frame* loc_slit_map_frame=NULL;
00391 const char* name=NULL;
00392 int size=0;
00393 int i=0;
00394 int idx=0;
00395 float lambda_step=0.1;
00396
00397 double *vlambda = NULL, *vslit = NULL, *vorder = NULL;
00398 double *pos_x = NULL, *pos_y = NULL;
00399 double j=0;
00400 double k=0;
00401
00402 xsh_wavemap_list *wm_list = NULL ;
00403 cpl_frame* result=NULL;
00404 char wm_tag[256];
00405 char wm_name[256];
00406
00407
00408 XSH_ASSURE_NOT_NULL_MSG(mod_cfg_frame,"Null model cfg frame!");
00409 XSH_ASSURE_NOT_NULL_MSG(spectral_format_frame,"Null spectral format frame!");
00410 XSH_ASSURE_NOT_NULL_MSG( instrument,"Null instrument setting!");
00411
00412 if (pre_frame != NULL){
00413 check( pre = xsh_pre_load( pre_frame, instrument));
00414 }
00415
00416
00417
00418
00419 check(spec_list=xsh_spectralformat_list_load(spectral_format_frame,
00420 instrument));
00421 check( xsh_model_config_load_best( mod_cfg_frame, &model_config));
00422 check( xsh_model_binxy( &model_config, instrument->binx,
00423 instrument->biny));
00424
00425 name=cpl_frame_get_filename(wave_map_frame);
00426 wmap_ima=cpl_image_load(name,CPL_TYPE_FLOAT,0,0);
00427
00428 name=cpl_frame_get_filename(slit_map_frame);
00429 smap_ima=cpl_image_load(name,CPL_TYPE_FLOAT,0,0);
00430
00431
00432
00433
00434 xsh_free_image(&wmap_ima);
00435
00436 sol_min_slit = cpl_image_get_min(smap_ima);
00437 sol_max_slit = cpl_image_get_max(smap_ima);
00438 xsh_free_image(&smap_ima);
00439
00440
00441 slit_size = (int)((sol_max_slit-sol_min_slit)/slit_step)+1;
00442
00443 size=0;
00444 for( i=0; i< spec_list->size; i++){
00445 float lambda_min, lambda_max;
00446 int lambda_size;
00447
00448 sol_min_order = (sol_min_order > spec_list->list[i].absorder) ? spec_list->list[i].absorder : sol_min_order;
00449 sol_max_order = (sol_max_order < spec_list->list[i].absorder) ? spec_list->list[i].absorder : sol_max_order;
00450
00451 lambda_min = spec_list->list[i].lambda_min_full;
00452 lambda_max = spec_list->list[i].lambda_max_full;
00453 XSH_ASSURE_NOT_ILLEGAL_MSG(lambda_max>= lambda_min,
00454 "lambda_max< lambda_min!! Check input spectralformat table, columns WLMAXFUL and WLMINFUL");
00455 lambda_size = (int)((lambda_max-lambda_min)/lambda_step)+1;
00456 size += lambda_size*slit_size;
00457 }
00458 xsh_msg_dbg_medium( "size %d", size ) ;
00459
00460
00461
00462
00463 idx = 0;
00464 XSH_MALLOC( vorder, double, size);
00465 XSH_MALLOC( vlambda, double, size);
00466 XSH_MALLOC( vslit, double, size);
00467 XSH_MALLOC( pos_x, double, size);
00468 XSH_MALLOC( pos_y, double, size);
00469
00470 for( i=0; i< spec_list->size; i++){
00471 double absorder;
00472 float lambda_min, lambda_max;
00473 int morder=0;
00474
00475 absorder= (double)spec_list->list[i].absorder;
00476 morder= spec_list->list[i].absorder;
00477 lambda_min = spec_list->list[i].lambda_min_full;
00478 lambda_max = spec_list->list[i].lambda_max_full;
00479 xsh_msg("order %f lambda %f-%f",absorder, lambda_min, lambda_max);
00480
00481 for(j=lambda_min; j <= lambda_max; j+=lambda_step){
00482 for( k=sol_min_slit; k <= sol_max_slit; k+=slit_step){
00483 double x, y;
00484 xsh_model_get_xy(&model_config,instrument,j,morder,k,&x,&y);
00485 vorder[idx] = absorder;
00486 vlambda[idx] = j;
00487 vslit[idx] = k;
00488 pos_x[idx] = x;
00489 pos_y[idx] = y;
00490 idx++;
00491 }
00492 }
00493 xsh_msg_dbg_medium( "i %d idx %d / %d", i, idx, size);
00494 }
00495
00496 if ( xsh_debug_level_get() >= XSH_DEBUG_LEVEL_MEDIUM){
00497 FILE *test = NULL;
00498 int itest;
00499
00500 test = fopen( "wavemap_grid.log", "w");
00501 for( itest=0; itest < size; itest++){
00502 fprintf(test, "%f %f\n", pos_x[itest], pos_y[itest]);
00503 }
00504 fclose(test);
00505 }
00506
00507
00508 check( wm_list = xsh_wavemap_list_create( instrument));
00509
00510 check( xsh_wavemap_list_compute_poly( vlambda, vslit, pos_x, pos_y, idx,
00511 vorder, dispsol_param, wm_list));
00512
00513 sprintf(wm_tag,"WAVE_MAP_POLY_%s",xsh_instrument_arm_tostring(instrument));
00514 sprintf(wm_name,"%s.fits",wm_tag);
00515 check( result = xsh_wavemap_list_save_poly( wm_list, order_tab_frame,
00516 pre, instrument, wm_tag,
00517 &disp_frame,&loc_slit_map_frame));
00518
00519 if(clean_tmp) {
00520 xsh_add_temporary_file(wm_name);
00521 }
00522
00523
00524
00525
00526
00527
00528
00529
00530
00531
00532 cleanup:
00533 xsh_free_image(&wmap_ima);
00534 xsh_free_image(&smap_ima);
00535 xsh_free_frame(&result);
00536 xsh_free_frame(&loc_slit_map_frame);
00537
00538 XSH_FREE( vorder);
00539 XSH_FREE( vlambda);
00540 XSH_FREE( vslit);
00541 XSH_FREE( pos_x);
00542 XSH_FREE( pos_y);
00543 xsh_pre_free( &pre);
00544 xsh_spectralformat_list_free( &spec_list);
00545 xsh_wavemap_list_free( &wm_list);
00546 return disp_frame;
00547 }
00548
00549
00550
00562
00563
00564 cpl_error_code
00565 xsh_wavemap_qc(cpl_frame* frm_map,const cpl_frame* frm_tab)
00566 {
00567 int sx=0;
00568
00569 double* cx=NULL;
00570 double* cy=NULL;
00571 int wx=0;
00572 int wy=0;
00573 cpl_image* ima=NULL;
00574 const char* name_tab=NULL;
00575 const char* name_map=NULL;
00576 char qc_wlen[40];
00577 double* pima=NULL;
00578 cpl_propertylist* plist=NULL;
00579 cpl_table* tab=NULL;
00580 cpl_table* ext=NULL;
00581
00582 double wlen=0;
00583 int ord_min=0;
00584 int ord_max=0;
00585 int i=0;
00586 int next=0;
00587
00588 XSH_ASSURE_NOT_NULL(frm_map);
00589 XSH_ASSURE_NOT_NULL(frm_tab);
00590 check(name_tab=cpl_frame_get_filename(frm_tab));
00591 check(tab=cpl_table_load(name_tab,2,0));
00592 check(ord_min=cpl_table_get_column_min(tab,"ABSORDER"));
00593 check(ord_max=cpl_table_get_column_max(tab,"ABSORDER"));
00594
00595
00596 name_map=cpl_frame_get_filename(frm_map);
00597 ima=cpl_image_load(name_map,CPL_TYPE_DOUBLE,0,0);
00598 pima=cpl_image_get_data_double(ima);
00599 sx=cpl_image_get_size_x(ima);
00600
00601 plist=cpl_propertylist_load(name_map,0);
00602 for(i=ord_min;i<=ord_max;i++) {
00603 next=cpl_table_and_selected_int(tab,"ABSORDER",CPL_EQUAL_TO,i);;
00604 ext=cpl_table_extract_selected(tab);
00605 cx=cpl_table_get_data_double(ext,"CENTER_X");
00606 cy=cpl_table_get_data_double(ext,"CENTER_Y");
00607 wx=cx[next/2];
00608 wy=cy[next/2];
00609 wlen=pima[wx+sx*wy];
00610 sprintf(qc_wlen,"%s%d",XSH_QC_WMAP_WAVEC,i);
00611 cpl_propertylist_append_double(plist,qc_wlen,wlen);
00612 xsh_free_table(&ext);
00613 cpl_table_select_all(tab);
00614
00615 }
00616
00617 check(xsh_update_pheader_in_image_multi(frm_map,plist));
00618
00619 cleanup:
00620 xsh_free_image(&ima);
00621 xsh_free_table(&tab);
00622 xsh_free_table(&ext);
00623 xsh_free_propertylist(&plist);
00624 return cpl_error_get_code();
00625 }
00626
00627
00628
00640
00641
00642 cpl_error_code
00643 xsh_wavetab_qc(cpl_frame* frm_tab, const int is_poly)
00644 {
00645
00646 cpl_table* ext=NULL;
00647 const char* name_tab=NULL;
00648 cpl_table* tab=NULL;
00649 cpl_table* tbl=NULL;
00650 int nsel=0;
00651 int ord_min=0;
00652 int ord_max=0;
00653 int i=0;
00654 int j=0;
00655 cpl_vector* loc_vec=NULL;
00656 double* pvec=NULL;
00657 double* pvec_all_ord=NULL;
00658 double* py=NULL;
00659 double ymin=0;
00660 double ymax=0;
00661 double ymin_all=FLT_MAX;
00662 double ymax_all=FLT_MIN;
00663
00664 double ymed=0;
00665 double yavg=0;
00666 char qc_line[40];
00667 cpl_propertylist* plist=NULL;
00668 cpl_vector* vec_all_ord=NULL;
00669 int ymin_all_ord=0;
00670 int ymax_all_ord=0;
00671
00672 XSH_ASSURE_NOT_NULL(frm_tab);
00673
00674 check(name_tab=cpl_frame_get_filename(frm_tab));
00675 check(tab=cpl_table_load(name_tab,1,0));
00676 check(ord_min=cpl_table_get_column_min(tab,"Order"));
00677 check(ord_max=cpl_table_get_column_max(tab,"Order"));
00678 check(plist=cpl_propertylist_load(name_tab,0));
00679
00680 check(nsel=cpl_table_and_selected_int(tab,"Slit_index",CPL_EQUAL_TO,4));
00681 check(ext=cpl_table_extract_selected(tab));
00682
00683 check(vec_all_ord=cpl_vector_new(ord_max-ord_min+1));
00684 check(pvec_all_ord=cpl_vector_get_data(vec_all_ord));
00685
00686 for(i=ord_min;i<=ord_max;i++) {
00687 check(nsel=cpl_table_and_selected_int(ext,"Order",CPL_EQUAL_TO,i));
00688 check(tbl=cpl_table_extract_selected(ext));
00689
00690 if(nsel>1) {
00691 check(loc_vec=cpl_vector_new(nsel-1));
00692 check(pvec=cpl_vector_get_data(loc_vec));
00693
00694 if(is_poly) {
00695 py=cpl_table_get_data_double(tbl,"Ypoly");
00696 } else {
00697 py=cpl_table_get_data_double(tbl,"Ythanneal");
00698 }
00699
00700 for(j=0;j<nsel-1;j++) {
00701 pvec[j]=fabs(py[j+1]-py[j]);
00702
00703 if(pvec[j]>ymax_all) {
00704 ymax_all=pvec[j];
00705 ymax_all_ord=i;
00706 }
00707 if(pvec[j]<ymin_all) {
00708 ymin_all=pvec[j];
00709 ymin_all_ord=i;
00710 }
00711
00712
00713 }
00714
00715 check(ymin=cpl_vector_get_min(loc_vec));
00716 check(ymax=cpl_vector_get_max(loc_vec));
00717 check(ymed=cpl_vector_get_median(loc_vec));
00718 check(yavg=cpl_vector_get_mean(loc_vec));
00719 check(pvec_all_ord[i-ord_min]=yavg);
00720
00721
00722
00723 sprintf(qc_line,"%s%d",XSH_QC_LINE_DIFMIN,i);
00724 check(cpl_propertylist_append_double(plist,qc_line,ymin));
00725 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMIN_C));
00726
00727 sprintf(qc_line,"%s%d",XSH_QC_LINE_DIFMAX,i);
00728 check(cpl_propertylist_append_double(plist,qc_line,ymax));
00729 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMAX_C));
00730
00731 sprintf(qc_line,"%s%d",XSH_QC_LINE_DIFMED,i);
00732 check(cpl_propertylist_append_double(plist,qc_line,ymed));
00733 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMED_C));
00734
00735 sprintf(qc_line,"%s%d",XSH_QC_LINE_DIFAVG,i);
00736 check(cpl_propertylist_append_double(plist,qc_line,yavg));
00737 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFAVG_C));
00738
00739 xsh_free_vector(&loc_vec);
00740 } else {
00741 xsh_msg_warning("Too few values of 'Slit_index=4' for Order=%d",i);
00742 xsh_msg_warning("No %s (and similar) QC parameters can be generated for order %d",XSH_QC_LINE_DIFMIN,i);
00743 }
00744 check(cpl_table_select_all(ext));
00745 xsh_free_table(&tbl);
00746 }
00747 check(yavg=cpl_vector_get_mean(vec_all_ord));
00748 check(ymax=cpl_vector_get_max(vec_all_ord));
00749 check(ymed=cpl_vector_get_median(vec_all_ord));
00750 check(ymin=cpl_vector_get_min(vec_all_ord));
00751
00752
00753
00754 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMIN);
00755 check(cpl_propertylist_append_double(plist,qc_line,ymin_all));
00756 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMIN_C));
00757
00758 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMIN_ORD);
00759 check(cpl_propertylist_append_int(plist,qc_line,ymin_all_ord));
00760 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMIN_C));
00761
00762 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMAX);
00763 check(cpl_propertylist_append_double(plist,qc_line,ymax_all));
00764 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMAX_C));
00765
00766 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMAX_ORD);
00767 check(cpl_propertylist_append_int(plist,qc_line,ymax_all_ord));
00768 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMAX_C));
00769
00770 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMED);
00771 check(cpl_propertylist_append_double(plist,qc_line,ymed));
00772 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMED_C));
00773
00774 sprintf(qc_line,"%s",XSH_QC_LINE_DIFAVG);
00775 check(cpl_propertylist_append_double(plist,qc_line,yavg));
00776 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFAVG_C));
00777
00778 check(cpl_table_select_all(tab));
00779 check(cpl_table_save(tab, plist, NULL,name_tab, CPL_IO_DEFAULT));
00780
00781 cleanup:
00782 xsh_free_table(&tab);
00783 xsh_free_table(&ext);
00784 xsh_free_table(&tbl);
00785 xsh_free_vector(&loc_vec);
00786 xsh_free_vector(&vec_all_ord);
00787 xsh_free_propertylist(&plist);
00788
00789
00790 return cpl_error_get_code();
00791
00792 }
00793