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[80];
00085 char slitmap_tag[80];
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 float sol_min_lambda=0, sol_max_lambda=0;
00213 float sol_min_slit, sol_max_slit;
00214 int sol_min_order, sol_max_order;
00215 int i, idx, size, slit_size;
00216 double j, k;
00217 double *vlambda = NULL, *vslit = NULL, *vorder = NULL;
00218 double *pos_x = NULL, *pos_y = NULL;
00219
00220 xsh_wavemap_list *wm_list = NULL ;
00221 int binx, biny ;
00222 char wm_name[80];
00223
00224
00225 XSH_ASSURE_NOT_NULL( wave_tab_2d_frame);
00226 XSH_ASSURE_NOT_NULL( order_tab_frame);
00227 XSH_ASSURE_NOT_NULL( spectral_format_frame);
00228 XSH_ASSURE_NOT_NULL( dispsol_par);
00229 XSH_ASSURE_NOT_NULL( instrument);
00230 XSH_ASSURE_NOT_NULL( dispersol_frame);
00231
00232
00233 if (pre_frame != NULL){
00234 check( pre = xsh_pre_load( pre_frame, instrument));
00235 }
00236
00237 check( binx = xsh_instrument_get_binx( instrument )) ;
00238 check( biny = xsh_instrument_get_biny( instrument )) ;
00239
00240 check( spec_list = xsh_spectralformat_list_load( spectral_format_frame,
00241 instrument));
00242 check( wave_tab_2d = xsh_wavesol_load( wave_tab_2d_frame,
00243 instrument ));
00244
00245 check( xsh_wavesol_set_bin_x( wave_tab_2d, binx ) ) ;
00246 check( xsh_wavesol_set_bin_y( wave_tab_2d, binx ) ) ;
00247
00248 sol_min_lambda = wave_tab_2d->min_lambda;
00249 sol_max_lambda = wave_tab_2d->max_lambda;
00250 sol_min_order = wave_tab_2d->min_order;
00251 sol_max_order = wave_tab_2d->max_order;
00252 sol_min_slit = wave_tab_2d->min_slit;
00253 sol_max_slit = wave_tab_2d->max_slit;
00254 xsh_msg("sol lambda %f,%f order %d,%d slit %f,%f",
00255 sol_min_lambda,sol_max_lambda, sol_min_order, sol_max_order,
00256 sol_min_slit, sol_max_slit);
00257
00258
00259 slit_size = (int)((sol_max_slit-sol_min_slit)/slit_step)+1;
00260
00261 size=0;
00262 for( i=0; i< spec_list->size; i++){
00263 float lambda_min, lambda_max;
00264 int lambda_size;
00265
00266 lambda_min = spec_list->list[i].lambda_min_full;
00267 lambda_max = spec_list->list[i].lambda_max_full;
00268 XSH_ASSURE_NOT_ILLEGAL_MSG(lambda_max>= lambda_min,
00269 "lambda_max< lambda_min!! Check input spectralformat table, columns WLMAXFUL and WLMINFUL");
00270 lambda_size = (int)((lambda_max-lambda_min)/lambda_step)+1;
00271 size += lambda_size*slit_size;
00272 }
00273 xsh_msg_dbg_medium( "size %d", size ) ;
00274
00275 idx = 0;
00276 XSH_MALLOC( vorder, double, size);
00277 XSH_MALLOC( vlambda, double, size);
00278 XSH_MALLOC( vslit, double, size);
00279 XSH_MALLOC( pos_x, double, size);
00280 XSH_MALLOC( pos_y, double, size);
00281
00282 for( i=0; i< spec_list->size; i++){
00283 double absorder;
00284 float lambda_min, lambda_max;
00285
00286 absorder= (double)spec_list->list[i].absorder;
00287 lambda_min = spec_list->list[i].lambda_min_full;
00288 lambda_max = spec_list->list[i].lambda_max_full;
00289 xsh_msg("order %f lambda %f-%f",absorder, lambda_min, lambda_max);
00290
00291 for(j=lambda_min; j <= lambda_max; j+=lambda_step){
00292 for( k=sol_min_slit; k <= sol_max_slit; k+=slit_step){
00293 double x, y;
00294
00295 check( x = xsh_wavesol_eval_polx( wave_tab_2d, j, absorder, k));
00296 check( y = xsh_wavesol_eval_poly( wave_tab_2d, j, absorder, k));
00297 vorder[idx] = absorder;
00298 vlambda[idx] = j;
00299 vslit[idx] = k;
00300 pos_x[idx] = x;
00301 pos_y[idx] = y;
00302 idx++;
00303 }
00304 }
00305 xsh_msg_dbg_medium( "i %d idx %d / %d", i, idx, size);
00306 }
00307
00308 if ( xsh_debug_level_get() >= XSH_DEBUG_LEVEL_MEDIUM){
00309 FILE *test = NULL;
00310 int itest;
00311
00312 test = fopen( "wavemap_grid.log", "w");
00313 for( itest=0; itest < size; itest++){
00314 fprintf(test, "%f %f\n", pos_x[itest], pos_y[itest]);
00315 }
00316 fclose(test);
00317 }
00318
00319
00320 check( wm_list = xsh_wavemap_list_create( instrument));
00321 #if 0
00322 check( xsh_wavemap_list_compute( vlambda, pos_x, pos_y, size,
00323 vorder, wm_par, wm_list));
00324 #else
00325
00326
00327
00328
00329
00330
00331 check( xsh_wavemap_list_compute_poly( vlambda, vslit, pos_x, pos_y, idx,
00332 vorder, dispsol_par, wm_list));
00333 #endif
00334 sprintf(wm_name,"%s.fits",wm_tag);
00335 check( result = xsh_wavemap_list_save_poly( wm_list, order_tab_frame,
00336 pre, instrument, wm_tag, dispersol_frame, slitmap_frame));
00337 if ( pre != NULL){
00338
00339 check( cpl_frame_set_tag( result,wm_tag));
00340 check( cpl_frame_set_tag( *slitmap_frame,
00341 XSH_GET_TAG_FROM_ARM( XSH_SLIT_MAP_POLY, instrument)));
00342 }
00343 cleanup:
00344
00345 XSH_FREE( vorder);
00346 XSH_FREE( vlambda);
00347 XSH_FREE( vslit);
00348 XSH_FREE( pos_x);
00349 XSH_FREE( pos_y);
00350 xsh_pre_free( &pre);
00351 xsh_spectralformat_list_free( &spec_list);
00352 xsh_wavesol_free( &wave_tab_2d);
00353 xsh_wavemap_list_free( &wm_list);
00354 return result;
00355 }
00356
00357
00358
00359 cpl_frame*
00360 xsh_create_dispersol_physmod(cpl_frame *pre_frame,
00361 cpl_frame *order_tab_frame,
00362 cpl_frame* mod_cfg_frame,
00363 cpl_frame* wave_map_frame,
00364 cpl_frame* slit_map_frame,
00365 xsh_dispersol_param *dispsol_param,
00366 cpl_frame* spectral_format_frame,
00367 xsh_instrument* instrument,
00368 const int clean_tmp)
00369 {
00370
00371
00372 xsh_pre *pre = NULL;
00373 cpl_frame* disp_frame=NULL;
00374 xsh_spectralformat_list *spec_list = NULL;
00375 int binx=0, biny=0 ;
00376 xsh_xs_3 model_config;
00377 int slit_size=0;
00378 int slit_step=1.5;
00379
00380
00381 float sol_min_lambda=0, sol_max_lambda=0;
00382 float sol_min_slit=0, sol_max_slit=0;
00383 int sol_min_order=0, sol_max_order=0;
00384 cpl_image* wmap_ima=NULL;
00385 cpl_image* smap_ima=NULL;
00386 const char* name=NULL;
00387 int size=0;
00388 int i=0;
00389 int idx=0;
00390 float lambda_step=0.1;
00391
00392 double *vlambda = NULL, *vslit = NULL, *vorder = NULL;
00393 double *pos_x = NULL, *pos_y = NULL;
00394 double j=0;
00395 double k=0;
00396
00397 xsh_wavemap_list *wm_list = NULL ;
00398 cpl_frame* result=NULL;
00399 char wm_tag[80];
00400 char wm_name[80];
00401
00402
00403 XSH_ASSURE_NOT_NULL_MSG(mod_cfg_frame,"Null model cfg frame!");
00404 XSH_ASSURE_NOT_NULL_MSG(spectral_format_frame,"Null spectral format frame!");
00405 XSH_ASSURE_NOT_NULL_MSG( instrument,"Null instrument setting!");
00406
00407 if (pre_frame != NULL){
00408 check( pre = xsh_pre_load( pre_frame, instrument));
00409 }
00410
00411
00412 check( binx = xsh_instrument_get_binx( instrument )) ;
00413 check( biny = xsh_instrument_get_biny( instrument )) ;
00414 check(spec_list=xsh_spectralformat_list_load(spectral_format_frame,
00415 instrument));
00416 check( xsh_model_config_load_best( mod_cfg_frame, &model_config));
00417 check( xsh_model_binxy( &model_config, instrument->binx,
00418 instrument->biny));
00419
00420 name=cpl_frame_get_filename(wave_map_frame);
00421 wmap_ima=cpl_image_load(name,CPL_TYPE_FLOAT,0,0);
00422
00423 name=cpl_frame_get_filename(slit_map_frame);
00424 smap_ima=cpl_image_load(name,CPL_TYPE_FLOAT,0,0);
00425
00426
00427 sol_min_lambda = cpl_image_get_min(wmap_ima);
00428 sol_max_lambda = cpl_image_get_max(wmap_ima);
00429 xsh_free_image(&wmap_ima);
00430
00431 sol_min_slit = cpl_image_get_min(smap_ima);
00432 sol_max_slit = cpl_image_get_max(smap_ima);
00433 xsh_free_image(&smap_ima);
00434
00435
00436 slit_size = (int)((sol_max_slit-sol_min_slit)/slit_step)+1;
00437
00438 size=0;
00439 for( i=0; i< spec_list->size; i++){
00440 float lambda_min, lambda_max;
00441 int lambda_size;
00442
00443 sol_min_order = (sol_min_order > spec_list->list[i].absorder) ? spec_list->list[i].absorder : sol_min_order;
00444 sol_max_order = (sol_max_order < spec_list->list[i].absorder) ? spec_list->list[i].absorder : sol_max_order;
00445
00446 lambda_min = spec_list->list[i].lambda_min_full;
00447 lambda_max = spec_list->list[i].lambda_max_full;
00448 XSH_ASSURE_NOT_ILLEGAL_MSG(lambda_max>= lambda_min,
00449 "lambda_max< lambda_min!! Check input spectralformat table, columns WLMAXFUL and WLMINFUL");
00450 lambda_size = (int)((lambda_max-lambda_min)/lambda_step)+1;
00451 size += lambda_size*slit_size;
00452 }
00453 xsh_msg_dbg_medium( "size %d", size ) ;
00454
00455
00456
00457
00458 idx = 0;
00459 XSH_MALLOC( vorder, double, size);
00460 XSH_MALLOC( vlambda, double, size);
00461 XSH_MALLOC( vslit, double, size);
00462 XSH_MALLOC( pos_x, double, size);
00463 XSH_MALLOC( pos_y, double, size);
00464
00465 for( i=0; i< spec_list->size; i++){
00466 double absorder;
00467 float lambda_min, lambda_max;
00468 int morder=0;
00469
00470 absorder= (double)spec_list->list[i].absorder;
00471 morder= spec_list->list[i].absorder;
00472 lambda_min = spec_list->list[i].lambda_min_full;
00473 lambda_max = spec_list->list[i].lambda_max_full;
00474 xsh_msg("order %f lambda %f-%f",absorder, lambda_min, lambda_max);
00475
00476 for(j=lambda_min; j <= lambda_max; j+=lambda_step){
00477 for( k=sol_min_slit; k <= sol_max_slit; k+=slit_step){
00478 double x, y;
00479 xsh_model_get_xy(&model_config,instrument,j,morder,k,&x,&y);
00480 vorder[idx] = absorder;
00481 vlambda[idx] = j;
00482 vslit[idx] = k;
00483 pos_x[idx] = x;
00484 pos_y[idx] = y;
00485 idx++;
00486 }
00487 }
00488 xsh_msg_dbg_medium( "i %d idx %d / %d", i, idx, size);
00489 }
00490
00491 if ( xsh_debug_level_get() >= XSH_DEBUG_LEVEL_MEDIUM){
00492 FILE *test = NULL;
00493 int itest;
00494
00495 test = fopen( "wavemap_grid.log", "w");
00496 for( itest=0; itest < size; itest++){
00497 fprintf(test, "%f %f\n", pos_x[itest], pos_y[itest]);
00498 }
00499 fclose(test);
00500 }
00501
00502
00503 check( wm_list = xsh_wavemap_list_create( instrument));
00504
00505 check( xsh_wavemap_list_compute_poly( vlambda, vslit, pos_x, pos_y, idx,
00506 vorder, dispsol_param, wm_list));
00507
00508 sprintf(wm_tag,"WAVE_MAP_POLY_%s",xsh_instrument_arm_tostring(instrument));
00509 sprintf(wm_name,"%s.fits",wm_tag);
00510 check( result = xsh_wavemap_list_save_poly( wm_list, order_tab_frame,
00511 pre, instrument, wm_tag,
00512 &disp_frame, &slit_map_frame));
00513
00514 if(clean_tmp) {
00515 xsh_add_temporary_file(wm_name);
00516 }
00517
00518
00519
00520
00521
00522
00523
00524
00525
00526
00527 cleanup:
00528 xsh_free_image(&wmap_ima);
00529 xsh_free_image(&smap_ima);
00530 xsh_free_frame(&result);
00531
00532 XSH_FREE( vorder);
00533 XSH_FREE( vlambda);
00534 XSH_FREE( vslit);
00535 XSH_FREE( pos_x);
00536 XSH_FREE( pos_y);
00537 xsh_pre_free( &pre);
00538 xsh_spectralformat_list_free( &spec_list);
00539 xsh_wavemap_list_free( &wm_list);
00540 return disp_frame;
00541 }
00542
00543
00544
00556
00557
00558 cpl_error_code
00559 xsh_wavemap_qc(cpl_frame* frm_map,const cpl_frame* frm_tab)
00560 {
00561 int sx=0;
00562 int sy=0;
00563 double* cx=NULL;
00564 double* cy=NULL;
00565 int wx=0;
00566 int wy=0;
00567 cpl_image* ima=NULL;
00568 const char* name_tab=NULL;
00569 const char* name_map=NULL;
00570 char qc_wlen[40];
00571 double* pima=NULL;
00572 cpl_propertylist* plist=NULL;
00573 cpl_table* tab=NULL;
00574 cpl_table* ext=NULL;
00575
00576 double wlen=0;
00577 int ord_min=0;
00578 int ord_max=0;
00579 int i=0;
00580 int next=0;
00581
00582 XSH_ASSURE_NOT_NULL(frm_map);
00583 XSH_ASSURE_NOT_NULL(frm_tab);
00584 check(name_tab=cpl_frame_get_filename(frm_tab));
00585 check(tab=cpl_table_load(name_tab,2,0));
00586 check(ord_min=cpl_table_get_column_min(tab,"ABSORDER"));
00587 check(ord_max=cpl_table_get_column_max(tab,"ABSORDER"));
00588
00589
00590 name_map=cpl_frame_get_filename(frm_map);
00591 ima=cpl_image_load(name_map,CPL_TYPE_DOUBLE,0,0);
00592 pima=cpl_image_get_data_double(ima);
00593 sx=cpl_image_get_size_x(ima);
00594 sy=cpl_image_get_size_y(ima);
00595 plist=cpl_propertylist_load(name_map,0);
00596 for(i=ord_min;i<=ord_max;i++) {
00597 next=cpl_table_and_selected_int(tab,"ABSORDER",CPL_EQUAL_TO,i);;
00598 ext=cpl_table_extract_selected(tab);
00599 cx=cpl_table_get_data_double(ext,"CENTER_X");
00600 cy=cpl_table_get_data_double(ext,"CENTER_Y");
00601 wx=cx[next/2];
00602 wy=cy[next/2];
00603 wlen=pima[wx+sx*wy];
00604 sprintf(qc_wlen,"%s%d",XSH_QC_WMAP_WAVEC,i);
00605 cpl_propertylist_append_double(plist,qc_wlen,wlen);
00606 xsh_free_table(&ext);
00607 cpl_table_select_all(tab);
00608
00609 }
00610
00611 check(xsh_update_pheader_in_image_multi(frm_map,plist));
00612
00613 cleanup:
00614 xsh_free_image(&ima);
00615 xsh_free_table(&tab);
00616 xsh_free_table(&ext);
00617 xsh_free_propertylist(&plist);
00618 return cpl_error_get_code();
00619 }
00620
00621
00622
00634
00635
00636 cpl_error_code
00637 xsh_wavetab_qc(cpl_frame* frm_tab, const int is_poly)
00638 {
00639
00640 cpl_table* ext=NULL;
00641 const char* name_tab=NULL;
00642 cpl_table* tab=NULL;
00643 cpl_table* tbl=NULL;
00644 int nsel=0;
00645 int ord_min=0;
00646 int ord_max=0;
00647 int i=0;
00648 int j=0;
00649 cpl_vector* loc_vec=NULL;
00650 double* pvec=NULL;
00651 double* pvec_all_ord=NULL;
00652 double* py=NULL;
00653 double ymin=0;
00654 double ymax=0;
00655 double ymin_all=FLT_MAX;
00656 double ymax_all=FLT_MIN;
00657
00658 double ymed=0;
00659 double yavg=0;
00660 char qc_line[40];
00661 cpl_propertylist* plist=NULL;
00662 cpl_vector* vec_all_ord=NULL;
00663 int ymin_all_ord=0;
00664 int ymax_all_ord=0;
00665
00666 XSH_ASSURE_NOT_NULL(frm_tab);
00667
00668 check(name_tab=cpl_frame_get_filename(frm_tab));
00669 check(tab=cpl_table_load(name_tab,1,0));
00670 check(ord_min=cpl_table_get_column_min(tab,"Order"));
00671 check(ord_max=cpl_table_get_column_max(tab,"Order"));
00672 check(plist=cpl_propertylist_load(name_tab,0));
00673
00674 check(nsel=cpl_table_and_selected_int(tab,"Slit_index",CPL_EQUAL_TO,4));
00675 check(ext=cpl_table_extract_selected(tab));
00676
00677 check(vec_all_ord=cpl_vector_new(ord_max-ord_min+1));
00678 check(pvec_all_ord=cpl_vector_get_data(vec_all_ord));
00679
00680 for(i=ord_min;i<=ord_max;i++) {
00681 check(nsel=cpl_table_and_selected_int(ext,"Order",CPL_EQUAL_TO,i));
00682 check(tbl=cpl_table_extract_selected(ext));
00683
00684 if(nsel>1) {
00685 check(loc_vec=cpl_vector_new(nsel-1));
00686 check(pvec=cpl_vector_get_data(loc_vec));
00687
00688 if(is_poly) {
00689 py=cpl_table_get_data_double(tbl,"Ypoly");
00690 } else {
00691 py=cpl_table_get_data_double(tbl,"Ythanneal");
00692 }
00693
00694 for(j=0;j<nsel-1;j++) {
00695 pvec[j]=fabs(py[j+1]-py[j]);
00696
00697 if(pvec[j]>ymax_all) {
00698 ymax_all=pvec[j];
00699 ymax_all_ord=i;
00700 }
00701 if(pvec[j]<ymin_all) {
00702 ymin_all=pvec[j];
00703 ymin_all_ord=i;
00704 }
00705
00706
00707 }
00708
00709 check(ymin=cpl_vector_get_min(loc_vec));
00710 check(ymax=cpl_vector_get_max(loc_vec));
00711 check(ymed=cpl_vector_get_median(loc_vec));
00712 check(yavg=cpl_vector_get_mean(loc_vec));
00713 check(pvec_all_ord[i-ord_min]=yavg);
00714
00715
00716
00717 sprintf(qc_line,"%s%d",XSH_QC_LINE_DIFMIN,i);
00718 check(cpl_propertylist_append_double(plist,qc_line,ymin));
00719 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMIN_C));
00720
00721 sprintf(qc_line,"%s%d",XSH_QC_LINE_DIFMAX,i);
00722 check(cpl_propertylist_append_double(plist,qc_line,ymax));
00723 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMAX_C));
00724
00725 sprintf(qc_line,"%s%d",XSH_QC_LINE_DIFMED,i);
00726 check(cpl_propertylist_append_double(plist,qc_line,ymed));
00727 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMED_C));
00728
00729 sprintf(qc_line,"%s%d",XSH_QC_LINE_DIFAVG,i);
00730 check(cpl_propertylist_append_double(plist,qc_line,yavg));
00731 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFAVG_C));
00732
00733 xsh_free_vector(&loc_vec);
00734 } else {
00735 xsh_msg_warning("Too few values of 'Slit_index=4' for Order=%d",i);
00736 xsh_msg_warning("No %s (and similar) QC parameters can be generated for order %d",XSH_QC_LINE_DIFMIN,i);
00737 }
00738 check(cpl_table_select_all(ext));
00739 xsh_free_table(&tbl);
00740 }
00741 check(yavg=cpl_vector_get_mean(vec_all_ord));
00742 check(ymax=cpl_vector_get_max(vec_all_ord));
00743 check(ymed=cpl_vector_get_median(vec_all_ord));
00744 check(ymin=cpl_vector_get_min(vec_all_ord));
00745
00746
00747
00748 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMIN);
00749 check(cpl_propertylist_append_double(plist,qc_line,ymin_all));
00750 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMIN_C));
00751
00752 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMIN_ORD);
00753 check(cpl_propertylist_append_int(plist,qc_line,ymin_all_ord));
00754 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMIN_C));
00755
00756 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMAX);
00757 check(cpl_propertylist_append_double(plist,qc_line,ymax_all));
00758 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMAX_C));
00759
00760 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMAX_ORD);
00761 check(cpl_propertylist_append_int(plist,qc_line,ymax_all_ord));
00762 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMAX_C));
00763
00764 sprintf(qc_line,"%s",XSH_QC_LINE_DIFMED);
00765 check(cpl_propertylist_append_double(plist,qc_line,ymed));
00766 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFMED_C));
00767
00768 sprintf(qc_line,"%s",XSH_QC_LINE_DIFAVG);
00769 check(cpl_propertylist_append_double(plist,qc_line,yavg));
00770 check(cpl_propertylist_set_comment(plist,qc_line,XSH_QC_LINE_DIFAVG_C));
00771
00772 check(cpl_table_select_all(tab));
00773 check(cpl_table_save(tab, plist, NULL,name_tab, CPL_IO_DEFAULT));
00774
00775 cleanup:
00776 xsh_free_table(&tab);
00777 xsh_free_table(&ext);
00778 xsh_free_table(&tbl);
00779 xsh_free_vector(&loc_vec);
00780 xsh_free_vector(&vec_all_ord);
00781 xsh_free_propertylist(&plist);
00782
00783
00784 return cpl_error_get_code();
00785
00786 }
00787