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
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037 #ifdef HAVE_CONFIG_H
00038 # include <config.h>
00039 #endif
00040
00041
00042
00043
00044 #include "sinfo_new_lamp_flats.h"
00045 #include "sinfo_flat_ini_by_cpl.h"
00046 #include "sinfo_pro_save.h"
00047 #include "sinfo_pro_types.h"
00048 #include "sinfo_functions.h"
00049 #include "sinfo_new_cube_ops.h"
00050 #include "sinfo_error.h"
00051 #include "sinfo_utils_wrappers.h"
00052 #include "sinfo_image_ops.h"
00053 #include "sinfo_utilities.h"
00054 #include "sinfo_globals.h"
00055
00056
00057
00058
00059 static int
00060 new_qc_get_cnt(cpl_frameset* on_set, cpl_frameset* of_set, flat_config* cfg);
00061 static int
00062 new_lamp_flats_det_ncounts(cpl_frameset* raw, flat_config* config);
00063
00064 static struct {
00065
00066 double avg_on;
00067 double std_on;
00068 double avg_of;
00069 double std_of;
00070 double avg_di;
00071 double std_di;
00072 double nsat_on;
00073 double noise_on;
00074 double noise_of;
00075 double flux_on;
00076 double nsat;
00077
00078 } qc_lampflat;
00079
00080
00081
00082
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112 int
00113 sinfo_new_lamp_flats (const char* plugin_id,
00114 cpl_parameterlist* config, cpl_frameset* sof)
00115 {
00116 flat_config * cfg =NULL;
00117 cpl_imagelist* list_object=NULL;
00118 cpl_imagelist* list_dither_object=NULL ;
00119 cpl_imagelist* list_sky=NULL ;
00120 cpl_imagelist* list_dither_sky=NULL;
00121 cpl_image ** im=NULL ;
00122 cpl_image * norm_dith =NULL;
00123 cpl_image * im_obj =NULL;
00124 cpl_image * int_im =NULL;
00125 cpl_image * int_im_dith =NULL;
00126 cpl_image * im_sky =NULL;
00127 cpl_image * im_dither =NULL;
00128 cpl_image * im_obj_sub =NULL;
00129 cpl_image * im_dither_sub =NULL;
00130 cpl_image * im_dither_sky =NULL;
00131 cpl_image ** imMed=NULL ;
00132 cpl_image * colImage =NULL;
00133 cpl_image * mask_im =NULL;
00134 cpl_image * norm =NULL;
00135 cpl_image * compImage =NULL;
00136 cpl_image * maskImage =NULL;
00137 cpl_image * threshIm =NULL;
00138
00139 char name[MAX_NAME_SIZE];
00140
00141 int typ;
00142 Stats * stats =NULL;
00143 int i = 0;
00144 int* n=NULL;
00145 int nob =0;
00146 int nsky = 0;
00147 int nobjdith = 0;
00148 int nskydith = 0;
00149 int n_badpixels =0;
00150 int pos =0;
00151
00152 float** slit_edges=NULL;
00153 float local_clean_mean =0.;
00154 float clean_stdev =0.;
00155 float mean_factor =0.;
00156 float val_x=0;
00157 float val_y=0;
00158
00159 char outNameDither[MAX_NAME_SIZE];
00160 char name_list[MAX_NAME_SIZE];
00161 char tbl_slitpos_name[MAX_NAME_SIZE];
00162
00163 cpl_table* tbl_slitpos=NULL;
00164 int* status=NULL;
00165
00166 int n_im_med=0;
00167 cpl_frameset* raw=NULL;
00168
00169 cpl_table* qclog_tbl=NULL;
00170 int naxis1=0;
00171 int naxis2=0;
00172 double fpn_stdev1=0;
00173 double fpn_stdev2=0;
00174
00175 int no=0;
00176 float lo_cut=0;
00177 float hi_cut=0;
00178
00179
00180
00181
00182
00183
00184
00185
00186 cknull_nomsg(raw=cpl_frameset_new());
00187
00188 cknull(cfg = sinfo_parse_cpl_input_flat(config,sof,&raw),
00189 "could not parse cpl input!");
00190
00191 if (cfg->interpolInd == 1) {
00192 if(sinfo_is_fits_file(cfg->mask) != 1) {
00193 sinfo_msg_error("Input file %s is not FITS",cfg->mask);
00194 goto cleanup;
00195 }
00196 if (sinfo_is_fits_file(cfg->slitposList) != 1) {
00197 sinfo_msg_error("Input file %s is not FITS",cfg->slitposList);
00198 goto cleanup;
00199 }
00200 }
00201
00202
00203
00204
00205
00206
00207
00208 sinfo_msg("Takes clean mean of several images");
00209
00210 cknull(list_object = cpl_imagelist_new (),"could not allocate memory");
00211
00212 if (cfg->contains_dither == 1) {
00213 cknull(list_dither_object=cpl_imagelist_new(),"could not allocate memory");
00214 }
00215
00216 if (cfg->contains_sky == 1) {
00217 cknull(list_sky=cpl_imagelist_new(),"could not allocate memory");
00218 }
00219
00220 if (cfg->contains_dither == 1 && cfg->nditheroff > 0) {
00221 cknull(list_dither_sky=cpl_imagelist_new(),"could not allocate memory");
00222 }
00223
00224 if (cfg->contains_dither == 0 && cfg->nditheroff > 0){
00225 sinfo_msg_error("please use non-dithered off-frames, remove the 2!");
00226 goto cleanup;
00227 }
00228
00229
00230
00231 im = (cpl_image**) cpl_calloc (cfg -> nframes, sizeof(cpl_image*));
00232
00233 for (i=0; i< cfg->nframes; i++) {
00234 strcpy(name,cfg->framelist[i]);
00235 if(sinfo_is_fits_file(name) != 1) {
00236 sinfo_msg_error("PP Input file %s %d is not FITS",name,i);
00237 goto cleanup;
00238 }
00239 im[i]=cpl_image_load(name,CPL_TYPE_FLOAT,0,0);
00240
00241 }
00242
00243 for (i=0; i< cfg->nframes; i++) {
00244 typ = cfg->frametype[i];
00245 pos = cfg->frameposition[i];
00246 cknull(im[i],"could not load image %d",i);
00247 if (pos == 2) {
00248 if (typ == 1) {
00249 cpl_imagelist_set( list_object, cpl_image_duplicate(im[i]), nob );
00250 nob = nob + 1;
00251 } else {
00252 cpl_imagelist_set( list_sky, cpl_image_duplicate(im[i]), nsky );
00253 nsky = nsky + 1 ;
00254 }
00255 } else {
00256 if (typ == 1) {
00257 cpl_imagelist_set(list_dither_object,
00258 cpl_image_duplicate(im[i]), nobjdith );
00259 nobjdith = nobjdith + 1;
00260 } else {
00261 cpl_imagelist_set( list_dither_sky,
00262 cpl_image_duplicate(im[i]), nskydith );
00263 nskydith = nskydith + 1 ;
00264 }
00265 }
00266 }
00267
00268
00269 if (nob != cfg->nobj || cfg->noff != nsky ||
00270 nobjdith != cfg->nditherobj || nskydith != cfg->nditheroff) {
00271 sinfo_msg_error("something is wrong with the number of "
00272 "the different types of frames");
00273 goto cleanup;
00274 }
00275
00276
00277 sinfo_msg("Creates and fills cubes with the different image lists");
00278 cknull(list_object,"could not create data cube!");
00279
00280 if (cfg->contains_dither == 1) {
00281 cknull(list_dither_object,"could not create data cube!");
00282 }
00283 if (cfg->contains_sky == 1 && nsky > 0) {
00284 cknull(list_sky,"could not create data cube!");
00285 }
00286
00287 if (cfg->contains_dither == 1 && nskydith > 0) {
00288 cknull(list_dither_sky,"could not create data cube!");
00289 }
00290
00291
00292
00293 sinfo_msg("Takes the average of the different cubes");
00294 if (cfg->loReject*cfg->nobj < 1. && cfg->hiReject *cfg->nobj < 1.) {
00295 cknull(im_obj = sinfo_new_average_cube_to_image(list_object ),
00296 "sinfo_averageCubeToImage failed" );
00297 } else {
00298
00299 no=cpl_imagelist_get_size(list_object);
00300 lo_cut=(floor)(cfg->loReject*no+0.5);
00301 hi_cut=(floor)(cfg->hiReject*no+0.5);
00302 cknull(im_obj=cpl_imagelist_collapse_minmax_create(list_object,
00303 lo_cut,
00304 hi_cut),
00305 "sinfo_average_with_rejection failed" );
00306 }
00307 sinfo_free_imagelist(&list_object);
00308
00309 if (cfg->contains_sky == 1) {
00310 if (cfg->loReject * nsky < 1. && cfg->hiReject * nsky < 1.) {
00311 cknull(im_sky = sinfo_new_average_cube_to_image(list_sky ),
00312 "sinfo_new_average_cube_to_image failed" );
00313 } else {
00314
00315 no=cpl_imagelist_get_size(list_sky);
00316 lo_cut=(floor)(cfg->loReject*no+0.5);
00317 hi_cut=(floor)(cfg->hiReject*no+0.5);
00318 cknull(im_sky=cpl_imagelist_collapse_minmax_create(list_sky,lo_cut,hi_cut),
00319 "sinfo_average_with_rejection failed" );
00320 }
00321 sinfo_free_imagelist(&list_sky);
00322 }
00323
00324 if (cfg->contains_dither == 1) {
00325 if (cfg->loReject*nobjdith < 1. && cfg->hiReject * nobjdith < 1.) {
00326 cknull(im_dither = sinfo_new_average_cube_to_image(list_dither_object ),
00327 "sinfo_new_average_cube_to_image failed" );
00328 } else {
00329
00330
00331 no=cpl_imagelist_get_size(list_dither_object);
00332 lo_cut=(floor)(cfg->loReject*no+0.5);
00333 hi_cut=(floor)(cfg->hiReject*no+0.5);
00334 cknull(im_dither=cpl_imagelist_collapse_minmax_create(list_dither_object,
00335 lo_cut,hi_cut),
00336 "sinfo_average_with_rejection failed" );
00337 }
00338 sinfo_free_imagelist(&list_dither_object);
00339 }
00340
00341 if (cfg->contains_dither == 1 && nskydith > 0 ) {
00342 if (cfg->loReject*nskydith < 1. && cfg->hiReject*nskydith < 1.) {
00343 cknull(im_dither_sky = sinfo_new_average_cube_to_image(list_dither_sky ),
00344 "sinfo_new_average_cube_to_image failed" );
00345 } else {
00346 no=cpl_imagelist_get_size(list_dither_sky);
00347 lo_cut=(floor)(cfg->loReject*no+0.5);
00348 hi_cut=(floor)(cfg->hiReject*no+0.5);
00349 cknull(im_dither_sky=cpl_imagelist_collapse_minmax_create(list_dither_sky,
00350 lo_cut,hi_cut),
00351 "new_average_with_rejection failed" );
00352 }
00353 sinfo_free_imagelist(&list_dither_sky);
00354 }
00355
00356
00357
00358
00359
00360
00361
00362
00363
00364 sinfo_msg("Subtracts the resulting off-frame (sky) from the on-frame");
00365 if (cfg->contains_sky == 1) {
00366 cknull(im_obj_sub = cpl_image_subtract_create(im_obj, im_sky),
00367 "could not sinfo_sub_image");
00368 sinfo_free_image(&im_obj);
00369 if (((cfg->contains_dither == 1) && (nskydith > 0)) ||
00370 (cfg->contains_dither == 0)) {
00371 sinfo_free_image(&im_sky);
00372 }
00373 im_obj = im_obj_sub;
00374 }
00375
00376 if (cfg->contains_dither == 1 && nskydith > 0) {
00377 cknull(im_dither_sub=cpl_image_subtract_create(im_dither, im_dither_sky),
00378 "could not sinfo_sub_image");
00379 sinfo_free_image(&im_dither);
00380 sinfo_free_image(&im_dither_sky);
00381 im_dither = im_dither_sub;
00382 } else if (cfg->contains_dither == 1 &&
00383 nskydith == 0 &&
00384 cfg->contains_sky == 1) {
00385 cknull(im_dither_sub = cpl_image_subtract_create(im_dither, im_sky),
00386 "could not sinfo_sub_image");
00387 sinfo_free_image(&im_dither);
00388 sinfo_free_image(&im_sky);
00389 im_dither = im_dither_sub;
00390 }
00391
00392
00393
00394
00395
00396
00397
00398
00399 sinfo_msg("Generating a static bad pixel mask");
00400 n_im_med = cfg->iterations+1;
00401
00402 imMed=(cpl_image**) cpl_calloc(n_im_med, sizeof(cpl_image*));
00403
00404 if (cfg->badInd == 1) {
00405 sinfo_msg("removes the intensity tilt from every column and");
00406 sinfo_msg("computes the standard deviation on a rectangular zone");
00407
00408
00409 cknull(colImage = sinfo_new_col_tilt( im_obj, cfg->sigmaFactor ),
00410 "sinfo_colTilt failed" );
00411
00412 cknull(stats = sinfo_new_image_stats_on_rectangle(colImage,
00413 cfg->badLoReject,
00414 cfg->badHiReject,
00415 cfg->llx,
00416 cfg->lly,
00417 cfg->urx,
00418 cfg->ury),
00419 "sinfo_get_image_stats_on_vig failed\n");
00420
00421 local_clean_mean = stats->cleanmean;
00422 clean_stdev = stats->cleanstdev;
00423
00424
00425
00426 if (cfg->threshInd == 1) {
00427 cknull(threshIm = sinfo_new_thresh_image(colImage,
00428 local_clean_mean-mean_factor*clean_stdev,
00429 local_clean_mean+mean_factor*clean_stdev),
00430 " sinfo_threshImage failed\n" );
00431 }
00432 if (cfg->threshInd == 0) {
00433 threshIm = colImage;
00434 }
00435
00436
00437
00438
00439
00440
00441
00442 cknull(imMed[0]= sinfo_new_median_image(threshIm,-cfg->factor*clean_stdev),
00443 " sinfo_medianImage failed" );
00444
00445
00446
00447
00448 for (i=1; i< cfg->iterations+1; i++) {
00449 cknull(imMed[i]=sinfo_new_median_image(imMed[i-1],
00450 -cfg->factor*clean_stdev),
00451 "sinfo_medianImage failed" );
00452 }
00453
00454
00455 cknull(compImage=sinfo_new_compare_images(threshIm,
00456 imMed[cfg->iterations],
00457 im_obj),
00458 "sinfo_compareImages failed" );
00459
00460
00461 n = (int*)cpl_calloc(1,sizeof(int));
00462 cknull(maskImage = sinfo_new_promote_image_to_mask( compImage, n ),
00463 "error in sinfo_promoteImageToMask" );
00464
00465
00466 n_badpixels = n[0];
00467 sinfo_msg("No of bad pixels: %d", n_badpixels);
00468
00469 cknull_nomsg(qclog_tbl = sinfo_qclog_init());
00470 ck0_nomsg(sinfo_qclog_add_int(qclog_tbl,"QC BP-MAP NBADPIX",n_badpixels,
00471 "No of bad pixels","%d"));
00472
00473 ck0(sinfo_pro_save_ima(maskImage,raw,sof,cfg->maskname,
00474 PRO_BP_MAP,qclog_tbl,plugin_id,config),
00475 "cannot save ima %s", cfg->maskname);
00476
00477
00478
00479 sinfo_free_table(&qclog_tbl);
00480
00481 sinfo_new_del_Stats(stats);
00482 sinfo_free_int(&n);
00483 sinfo_free_image(&threshIm);
00484 if (cfg->threshInd == 1) {
00485 sinfo_free_image(&colImage);
00486 }
00487 sinfo_free_image(&compImage);
00488 sinfo_free_image(&maskImage);
00489
00490 for (i=0; i< cfg->iterations+1; i++) {
00491 sinfo_free_image(&imMed[i]);
00492 }
00493
00494 }
00495
00496 cpl_free(imMed);
00497 imMed=NULL;
00498
00499
00500
00501
00502
00503
00504
00505 sinfo_msg("Creates a Master flat field");
00506 if (cfg->interpolInd == 1) {
00507 cknull(mask_im = cpl_image_load(cfg->mask,CPL_TYPE_FLOAT,0,0),
00508 "could not load static bad pixel mask" );
00509
00510
00511
00512
00513 slit_edges = (float **) cpl_calloc( 32, sizeof (float*) ) ;
00514 for ( i = 0 ; i < 32 ; i++ )
00515 {
00516 slit_edges[i] = (float *) cpl_calloc( 2, sizeof (float)) ;
00517 }
00518
00519 if(sinfo_is_fits_file(cfg->slitposList) !=1 ) {
00520 sinfo_msg_error("Input file %s is not FITS", cfg->slitposList);
00521 goto cleanup;
00522 }
00523 strcpy(tbl_slitpos_name,cfg->slitposList);
00524 check(tbl_slitpos = cpl_table_load(tbl_slitpos_name,1,0),
00525 "error loading tbl %s",tbl_slitpos_name);
00526
00527 for (i =0 ; i< 32; i++){
00528 val_x=cpl_table_get_double(tbl_slitpos,"pos1",i,status);
00529 val_y=cpl_table_get_double(tbl_slitpos,"pos2",i,status);
00530 slit_edges[i][0]=val_x;
00531 slit_edges[i][1]=val_y;
00532 }
00533 sinfo_free_table(&tbl_slitpos);
00534
00535 cknull(int_im = sinfo_interpol_source_image (im_obj, mask_im,
00536 cfg->maxRad, slit_edges),
00537 "could not carry out sinfo_interpolSourceImage" );
00538
00539 sinfo_free_image(&im_obj);
00540 cknull(norm = sinfo_new_normalize_to_central_pixel(int_im),
00541 "could not normalize flatfield" );
00542 sinfo_free_image(&int_im);
00543 im_obj = norm;
00544
00545 if (cfg->contains_dither == 1) {
00546 cknull(int_im_dith = sinfo_interpol_source_image(im_dither,
00547 mask_im,
00548 cfg->maxRad,
00549 slit_edges),
00550 "could not carry out sinfo_interpolSourceImage" );
00551
00552 cpl_image_delete(im_dither);
00553 cknull(norm_dith = sinfo_new_normalize_to_central_pixel(int_im_dith),
00554 "could not normalize flatfield" );
00555 sinfo_free_image(&int_im_dith);
00556 im_dither = norm_dith;
00557 }
00558
00559 for ( i = 0 ; i < 32 ; i++ )
00560 {
00561 cpl_free( slit_edges[i] );
00562 }
00563 cpl_free( slit_edges ) ;
00564 sinfo_free_image(&mask_im);
00565
00566 }
00567
00568 if (cfg->interpolInd != 1) {
00569 cknull(norm = sinfo_new_normalize_to_central_pixel(im_obj),
00570 "could not normalize flatfield" );
00571 sinfo_free_image(&im_obj);
00572 im_obj = norm;
00573
00574 if (cfg->contains_dither == 1) {
00575 cknull(norm_dith = sinfo_new_normalize_to_central_pixel(im_dither),
00576 "could not normalize flatfield" );
00577 sinfo_free_image(&im_dither);
00578 im_dither = norm_dith;
00579 }
00580 }
00581
00582 naxis1=cpl_image_get_size_x(im_obj);
00583 naxis2=cpl_image_get_size_y(im_obj);
00584
00585
00586 if(cfg->qc_fpn_xmin1 < 1) {
00587 sinfo_msg_error("qc_ron_xmin < 1");
00588 goto cleanup;
00589 }
00590
00591 if(cfg->qc_fpn_xmax1 > naxis1) {
00592 sinfo_msg_error("qc_ron_xmax < %d",naxis1);
00593 goto cleanup;
00594 }
00595
00596 if(cfg->qc_fpn_ymin1 < 1) {
00597 sinfo_msg_error("qc_ron_ymin < 1");
00598 goto cleanup;
00599 }
00600
00601 if(cfg->qc_fpn_ymax1 > naxis2) {
00602 sinfo_msg_error("qc_ron_ymax < %d",naxis2);
00603 goto cleanup;
00604 }
00605 fpn_stdev1 = cpl_image_get_stdev_window(im_obj,
00606 cfg->qc_fpn_xmin1,
00607 cfg->qc_fpn_ymin1,
00608 cfg->qc_fpn_xmax1,
00609 cfg->qc_fpn_ymax1);
00610
00611
00612 if(cfg->qc_fpn_xmin2 < 1) {
00613 sinfo_msg_error("qc_ron_xmin < %d",1);
00614 goto cleanup;
00615 }
00616
00617
00618 if(cfg->qc_fpn_xmax2 > naxis1) {
00619 sinfo_msg_error("qc_ron_xmax < %d",naxis1);
00620 goto cleanup;
00621 }
00622
00623 if(cfg->qc_fpn_ymin2 < 1) {
00624 sinfo_msg_error("qc_ron_ymin < 1");
00625 goto cleanup;
00626 }
00627
00628 if(cfg->qc_fpn_ymax2 > naxis2) {
00629 sinfo_msg_error("qc_ron_ymax < %d",naxis2);
00630 goto cleanup;
00631 }
00632 fpn_stdev2 = cpl_image_get_stdev_window(im_obj,
00633 cfg->qc_fpn_xmin2,
00634 cfg->qc_fpn_ymin2,
00635 cfg->qc_fpn_xmax2,
00636 cfg->qc_fpn_ymax2);
00637
00638
00639
00640 ck0(new_lamp_flats_det_ncounts(raw,cfg),"error computing number of counts");
00641 cknull_nomsg(qclog_tbl = sinfo_qclog_init());
00642 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC SPECFLAT NCNTSAVG",
00643 qc_lampflat.avg_di,"Average counts","%g"));
00644 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC SPECFLAT NCNTSSTD",
00645 qc_lampflat.std_di,"Stdev counts","%g"));
00646 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC SPECFLAT OFFFLUX",
00647 qc_lampflat.avg_of,
00648 "Average flux off frames","%g"));
00649
00650 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,
00651 "QC LFLAT FPN1",
00652 fpn_stdev1,
00653 "Fixed Pattern Noise of combined frames",
00654 "%f"));
00655
00656 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,
00657 "QC LFLAT FPN2",
00658 fpn_stdev2,
00659 "Fixed Pattern Noise of combined frames",
00660 "%f"));
00661
00662 ck0(sinfo_pro_save_ima(im_obj,raw,sof,cfg->outName,
00663 PRO_MASTER_FLAT_LAMP,qclog_tbl,plugin_id,config),
00664 "cannot save ima %s", cfg->outName);
00665
00666
00667 sinfo_free_table(&qclog_tbl);
00668 sinfo_free_image(&im_obj);
00669
00670
00671 if (cfg->contains_dither == 1) {
00672
00673 if (strstr(cfg->outName, ".fits" ) != NULL ) {
00674
00675 snprintf(name_list, MAX_NAME_SIZE-1,"%s%s", sinfo_new_get_rootname(cfg->outName),
00676 "_dither");
00677 strcpy(outNameDither,name_list);
00678 strcat(outNameDither,strstr(cfg->outName,".fits"));
00679
00680 } else {
00681 strcpy(outNameDither,cfg->outName);
00682 strcat(outNameDither,"_dither");
00683 }
00684
00685
00686 naxis1=cpl_image_get_size_x(im_dither);
00687 naxis2=cpl_image_get_size_y(im_dither);
00688
00689
00690 if(cfg->qc_fpn_xmin1 < 1) {
00691 sinfo_msg_error("qc_ron_xmin1 < 1");
00692 goto cleanup;
00693 }
00694
00695 if(cfg->qc_fpn_xmax1 > naxis1) {
00696 sinfo_msg_error("qc_ron_xmax1 < %d",naxis1);
00697 goto cleanup;
00698 }
00699
00700 if(cfg->qc_fpn_ymin1 < 1) {
00701 sinfo_msg_error("qc_ron_ymin1 < 1");
00702 goto cleanup;
00703 }
00704
00705 if(cfg->qc_fpn_ymax1 > naxis2) {
00706 sinfo_msg_error("qc_ron_ymax1 < %d",naxis2);
00707 goto cleanup;
00708 }
00709
00710
00711 fpn_stdev1 = cpl_image_get_stdev_window(im_dither,
00712 cfg->qc_fpn_xmin1,
00713 cfg->qc_fpn_ymin1,
00714 cfg->qc_fpn_xmax1,
00715 cfg->qc_fpn_ymax1);
00716
00717 if(cfg->qc_fpn_xmin2 < 1) {
00718 sinfo_msg_error("qc_ron_xmin2 < 1");
00719 goto cleanup;
00720 }
00721
00722 if(cfg->qc_fpn_xmax2 > naxis1) {
00723 sinfo_msg_error("qc_ron_xmax2 < %d",naxis1);
00724 goto cleanup;
00725 }
00726
00727 if(cfg->qc_fpn_ymin2 < 1) {
00728 sinfo_msg_error("qc_ron_ymin2 < 1");
00729 goto cleanup;
00730 }
00731
00732 if(cfg->qc_fpn_ymax2 > naxis2) {
00733 sinfo_msg_error("qc_ron_ymax2 < %d",naxis2);
00734 goto cleanup;
00735 }
00736
00737 fpn_stdev2 = cpl_image_get_stdev_window(im_dither,
00738 cfg->qc_fpn_xmin2,
00739 cfg->qc_fpn_ymin2,
00740 cfg->qc_fpn_xmax2,
00741 cfg->qc_fpn_ymax2);
00742
00743
00744 ck0(new_lamp_flats_det_ncounts(raw,cfg),"error computing ncounts");
00745 cknull_nomsg(qclog_tbl = sinfo_qclog_init());
00746 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC SPECFLAT NCNTSAVG",
00747 qc_lampflat.avg_di,"Average counts","%g"));
00748
00749 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC SPECFLAT NCNTSSTD",
00750 qc_lampflat.std_di,"Stdev counts","%g"));
00751
00752 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC SPECFLAT OFFFLUX",
00753 qc_lampflat.avg_of,"Average flux off frames","%g"));
00754
00755 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC LFLAT FPN1",fpn_stdev1,
00756 "Fixed Pattern Noise of combined frames","%f"));
00757
00758 ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC LFLAT FPN2",fpn_stdev2,
00759 "Fixed Pattern Noise of combined frames","%f"));
00760
00761
00762 ck0(sinfo_pro_save_ima(im_dither,raw,sof,outNameDither,
00763 PRO_MASTER_FLAT_LAMP,qclog_tbl,plugin_id,config),
00764 "cannot save ima %s", outNameDither);
00765
00766
00767 sinfo_free_table(&qclog_tbl);
00768 sinfo_free_image(&im_dither);
00769
00770 }
00771
00772
00773
00774 sinfo_free_image_array(&im,cfg->nframes);
00775 sinfo_free_frameset(&raw);
00776 sinfo_flat_free(&cfg);
00777 return 0;
00778
00779 cleanup:
00780
00781
00782
00783 if(slit_edges != NULL) {
00784 for ( i = 0 ; i < 32 ; i++ )
00785 {
00786 if(slit_edges[i] != NULL) {
00787 cpl_free( slit_edges[i] );
00788 }
00789 slit_edges[i]=NULL;
00790 }
00791 cpl_free( slit_edges ) ;
00792 }
00793 sinfo_free_image(&mask_im);
00794 sinfo_free_table(&qclog_tbl);
00795 if(stats) {
00796 sinfo_new_del_Stats(stats);
00797 }
00798 sinfo_free_int(&n);
00799 sinfo_free_image(&threshIm);
00800 sinfo_free_image(&maskImage);
00801 if(imMed != NULL) sinfo_free_image_array(&imMed,cfg->iterations);
00802 if(stats!= NULL) {
00803 sinfo_new_del_Stats(stats);
00804 stats=NULL;
00805 }
00806 sinfo_free_image(&compImage);
00807 sinfo_free_image(&colImage);
00808 sinfo_free_imagelist(&list_dither_object);
00809 sinfo_free_imagelist(&list_dither_sky);
00810 if(list_sky != NULL) {
00811 sinfo_free_imagelist(&list_sky);
00812 }
00813 sinfo_free_imagelist(&list_object);
00814 sinfo_free_image(&im_dither);
00815 sinfo_free_image(&im_dither_sky);
00816 sinfo_free_image(&im_obj);
00817 sinfo_free_image(&im_sky);
00818 if(im != NULL) sinfo_free_image_array(&im,cfg->nframes);
00819 sinfo_free_frameset(&raw);
00820 if(cfg != NULL) {
00821 sinfo_flat_free(&cfg);
00822 }
00823 return -1;
00824
00825 }
00826
00827 static int
00828 new_lamp_flats_det_ncounts(cpl_frameset* raw, flat_config* cfg)
00829 {
00830 int i=0;
00831 int j=0;
00832
00833 int nraw=0;
00834 int non=0;
00835 int noff=0;
00836
00837 double mjd_on=0;
00838 double mjd_of=0;
00839 double mjd_of_frm=0;
00840
00841 char filename[MAX_NAME_SIZE];
00842
00843 cpl_frame* frm=NULL;
00844 cpl_frame* frm_dup=NULL;
00845 cpl_frame* on_frm=NULL;
00846 cpl_frame* of_frm=NULL;
00847 cpl_frame* tmp_of_frm=NULL;
00848
00849
00850 cpl_frameset* on_set=NULL;
00851 cpl_frameset* of_set=NULL;
00852 cpl_frameset* wrk_set=NULL;
00853
00854 on_set=cpl_frameset_new();
00855 of_set=cpl_frameset_new();
00856
00857 nraw = cpl_frameset_get_size(raw);
00858
00859 for (i=0; i< nraw; i++) {
00860 frm = cpl_frameset_get_frame(raw,i);
00861 frm_dup = cpl_frame_duplicate(frm);
00862 if(sinfo_frame_is_on(frm) == 1) {
00863 cpl_frameset_insert(on_set,frm_dup);
00864 non++;
00865 } else {
00866 cpl_frameset_insert(of_set,frm_dup);
00867 noff++;
00868 }
00869 }
00870
00871
00872 if (non == noff) {
00873 new_qc_get_cnt(on_set,of_set,cfg);
00874
00875 } else if (non == 0) {
00876 sinfo_msg("non == 0");
00877 sinfo_msg_warning("QC SPECFLAT NCNTAVG/NCTNTSTD/OFFFLUX=0 ");
00878
00879 } else if ( noff == 0 ) {
00880 sinfo_msg("noff == 0");
00881 sinfo_msg_warning("QC SPECFLAT NCNTAVG/NCTNTSTD/OFFFLUX=0 ");
00882
00883 } else {
00884
00885 sinfo_msg_warning("non != noff, => QC SPECFLAT NCNTAVG/NCTNTSTD/OFFFLUX=0 ");
00886
00887 for (i=0;i<non;i++) {
00888 wrk_set=cpl_frameset_new();
00889 on_frm=cpl_frameset_get_frame(on_set,i);
00890 mjd_on=sinfo_get_mjd_obs(on_frm);
00891 of_frm=cpl_frameset_get_frame(of_set,0);
00892 mjd_of=sinfo_get_mjd_obs(of_frm);
00893 strcpy(filename,cpl_frame_get_filename(of_frm));
00894 for (j=1;j<noff;j++) {
00895 tmp_of_frm = cpl_frameset_get_frame(of_set,j);
00896 mjd_of_frm = sinfo_get_mjd_obs(tmp_of_frm);
00897
00898 if(1000.*(mjd_of_frm-mjd_on)*(mjd_of_frm-mjd_on) <
00899 1000.*(mjd_of- mjd_on)*(mjd_of- mjd_on) ) {
00900 mjd_of=mjd_of_frm;
00901 of_frm=cpl_frame_duplicate(tmp_of_frm);
00902 }
00903 }
00904 strcpy(filename,cpl_frame_get_filename(of_frm));
00905 frm_dup=cpl_frame_duplicate(of_frm);
00906 cpl_frameset_insert(wrk_set,frm_dup);
00907 strcpy(filename,cpl_frame_get_filename(of_frm));
00908 }
00909
00910 new_qc_get_cnt(on_set,wrk_set,cfg);
00911
00912 }
00913
00914 cpl_frameset_delete(wrk_set);
00915 cpl_frameset_delete(on_set);
00916 cpl_frameset_delete(of_set);
00917 return 0;
00918
00919
00920 }
00921
00922 static int
00923 new_qc_get_cnt(cpl_frameset* on_set, cpl_frameset* of_set, flat_config* cfg)
00924 {
00925
00926 int i=0;
00927 int nsat=0;
00928 int non=0;
00929 int nof=0;
00930 int nfr=0;
00931
00932 char name[MAX_NAME_SIZE];
00933 cpl_vector* vec_on=NULL;
00934 cpl_vector* vec_of=NULL;
00935 cpl_vector* vec_di=NULL;
00936 cpl_vector* vec_nsat=NULL;
00937 cpl_frame* on_frm=NULL;
00938 cpl_frame* of_frm=NULL;
00939
00940 cpl_image* dif_ima=NULL;
00941 cpl_image* on_ima=NULL;
00942 cpl_image* of_ima=NULL;
00943 cpl_image* tmp_ima=NULL;
00944
00945 double med=0;
00946 non = cpl_frameset_get_size(on_set);
00947 nof = cpl_frameset_get_size(of_set);
00948 nfr = (non <= nof) ? non : nof;
00949 vec_on = cpl_vector_new(nfr);
00950 vec_of = cpl_vector_new(nfr);
00951 vec_di = cpl_vector_new(nfr);
00952 vec_nsat = cpl_vector_new(nfr);
00953
00954
00955 for (i=0; i< nfr; i++) {
00956 on_frm = cpl_frameset_get_frame(on_set,i);
00957 strcpy(name,cpl_frame_get_filename(on_frm));
00958 on_ima = cpl_image_load(name,CPL_TYPE_FLOAT,0,0);
00959 med= cpl_image_get_median(on_ima);
00960 cpl_vector_set(vec_on,i,med);
00961
00962 tmp_ima = cpl_image_duplicate(on_ima);
00963 cpl_image_threshold(tmp_ima,CPL_PIXEL_MINVAL,
00964 cfg->qc_thresh_max,0,1);
00965 nsat=cpl_image_get_flux(tmp_ima);
00966 cpl_vector_set(vec_nsat,i,nsat);
00967
00968
00969 of_frm = cpl_frameset_get_frame(of_set,i);
00970 strcpy(name,cpl_frame_get_filename(of_frm));
00971 of_ima = cpl_image_load(name,CPL_TYPE_FLOAT,0,0);
00972 med= cpl_image_get_median(of_ima);
00973 cpl_vector_set(vec_of,i,med);
00974 dif_ima = cpl_image_subtract_create(on_ima,of_ima);
00975 med= cpl_image_get_median(dif_ima);
00976 cpl_vector_set(vec_di,i,med);
00977
00978 cpl_image_delete(on_ima);
00979 cpl_image_delete(of_ima);
00980 cpl_image_delete(dif_ima);
00981 cpl_image_delete(tmp_ima);
00982 }
00983 qc_lampflat.avg_on=cpl_vector_get_mean(vec_on);
00984 qc_lampflat.avg_of=cpl_vector_get_mean(vec_of);
00985 qc_lampflat.avg_di=cpl_vector_get_mean(vec_di);
00986 if(nfr > 1 ) {
00987 qc_lampflat.std_on=cpl_vector_get_stdev(vec_on);
00988 qc_lampflat.std_of=cpl_vector_get_stdev(vec_of);
00989 qc_lampflat.std_di=cpl_vector_get_stdev(vec_di);
00990 }
00991 qc_lampflat.nsat=cpl_vector_get_mean(vec_nsat);
00992 cpl_vector_delete(vec_on);
00993 cpl_vector_delete(vec_of);
00994 cpl_vector_delete(vec_di);
00995 cpl_vector_delete(vec_nsat);
00996
00997
00998
00999
01000
01001
01002
01003
01004
01005 return 0;
01006 }
01007