00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #include "twiflat.h"
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041 void change_plist_twiflat (cpl_propertylist* plist, char* out_name, int n_frames);
00042
00043 int twiflat (cpl_parameterlist* config,cpl_frameset* sof)
00044 {
00045 const char* _id = "twiflat";
00046 twiflat_config * cfg=NULL ;
00047
00048 OneImage ** list_twi ;
00049 OneImage ** list_dark ;
00050 OneImage ** list_on ;
00051 OneImage ** list_off ;
00052 OneImage ** im ;
00053
00054 OneImage * im_twi=NULL ;
00055 OneImage * im_dark=NULL ;
00056 OneImage * im_on=NULL ;
00057 OneImage * im_off=NULL ;
00058 OneImage * im_twi_sub=NULL ;
00059 OneImage * im_twi_bad=NULL ;
00060 OneImage * im_twi_tilt=NULL ;
00061 OneImage * im_on_sub=NULL ;
00062 OneImage * im_on_bad=NULL ;
00063 OneImage * im_on_tilt=NULL ;
00064
00065 OneImage * im_on_norm ;
00066
00067 OneImage * flat_field=NULL;
00068 OneImage * mask_im=NULL ;
00069 OneImage* row_twi=NULL;
00070 OneImage* row_on=NULL;
00071
00072 Vector* row_twi_vec=NULL;
00073 Vector* row_on_vec=NULL;
00074
00075 OneCube * cube_dark=NULL ;
00076 OneCube * cube_on=NULL ;
00077 OneCube * cube_off=NULL ;
00078 OneCube * cube_twi=NULL ;
00079
00080
00081
00082
00083 cpl_propertylist* cplist=NULL;
00084 char* name=NULL;
00085
00086
00087 int i = 0;
00088 int typ=0;
00089
00090 int ntwi=0;
00091 int ndark=0;
00092 int non=0;
00093 int noff=0;
00094
00095
00096 cpl_image* ff_img=NULL;
00097 cpl_frameset* raw=NULL;
00098
00099
00100
00101
00102
00103
00104
00105
00106 cfg = parse_cpl_input_twiflat(config,sof,&raw) ;
00107
00108 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00109 cpl_msg_error(_id,"error parsing cpl input");
00110 cpl_msg_error(_id,(char* ) cpl_error_get_message());
00111 return -1;
00112 }
00113 if (cfg == NULL)
00114 {
00115 cpl_msg_error (_id,"could not parse cpl input!\n") ;
00116 return -1 ;
00117 }
00118
00119
00120
00121 if(is_fits_file(cfg->mask) !=1 ){
00122 cpl_msg_error(_id,"Input file %s cfg->mask is not FITS",cfg->mask);
00123 return -1;
00124 }
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136
00137 list_twi = new_image_list (cfg->ntwi);
00138 if (list_twi == NULL) {
00139 cpl_msg_error(_id,"could not allocate memory\n");
00140 return -1;
00141 }
00142
00143 list_dark = new_image_list(cfg->ndark);
00144 if (list_dark == NULL) {
00145 cpl_msg_error(_id,"could not allocate memory\n");
00146 return -1;
00147 }
00148 list_on = new_image_list(cfg->non);
00149
00150 if (list_on == NULL) {
00151 cpl_msg_error(_id,"could not allocate memory\n");
00152 return -1;
00153 }
00154
00155 list_off = new_image_list(cfg->noff);
00156 if (list_off == NULL) {
00157 cpl_msg_error(_id,"could not allocate memory\n");
00158 return -1;
00159 }
00160
00161
00162 ntwi = 0;
00163 ndark = 0;
00164 non = 0;
00165 noff = 0;
00166
00167
00168
00169 im = (OneImage**) cpl_calloc(cfg -> nframes, sizeof(OneImage*));
00170 for (i=0; i< cfg->nframes; i++) {
00171 name = cfg->framelist[i];
00172 if(is_fits_file(name) != 1) {
00173 cpl_msg_info(_id,"Input file %s is not FITS",name);
00174 return -1;
00175 }
00176 im[i]=load_image( name );
00177
00178
00179
00180 }
00181
00182
00183 for (i=0; i< cfg->nframes; i++) {
00184 typ = intarray_get_value( cfg->frametype, i );
00185 if (im[i] == NULL) {
00186 cpl_msg_error(_id,"could not load image");
00187 return -1;
00188 }
00189 if (typ == 0) {
00190 insert_image( list_twi, im[i], ntwi );
00191
00192 ntwi = ntwi + 1;
00193 }
00194 else if (typ == 1) {
00195 insert_image( list_dark, im[i], ndark );
00196 ndark = ndark + 1;
00197 }
00198 else if (typ == 2) {
00199 insert_image( list_on, im[i], non );
00200 non = non + 1;
00201 }
00202 else if (typ == 3) {
00203 insert_image( list_off, im[i], noff );
00204 noff = noff + 1;
00205 }
00206 else {
00207 cpl_msg_error(_id,"wrong frame type\n");
00208 return -1;
00209 }
00210 }
00211 if ((cfg->non != non) ||
00212 (cfg->noff != noff) ||
00213 (cfg->ntwi != ntwi) ||
00214 (cfg->ndark != ndark)) {
00215 cpl_msg_error(_id,"something wrong with the number of the different types of frames");
00216 return -1;
00217 }
00218
00219
00220
00221 cube_twi = list_make_cube(list_twi, ntwi);
00222 if (cube_twi == NULL) {
00223 cpl_msg_error (_id,"could not create data cube!");
00224 return -1;
00225 }
00226
00227 destroy_list(list_twi);
00228 cube_dark = list_make_cube(list_dark, ndark);
00229
00230 if (cube_dark == NULL) {
00231 cpl_msg_error (_id,"could not create data cube!");
00232 return -1;
00233 }
00234
00235 destroy_list(list_dark);
00236 cube_on = list_make_cube(list_on, non);
00237 if (cube_on == NULL) {
00238 cpl_msg_error (_id,"could not create data cube!");
00239 return -1;
00240 }
00241
00242 destroy_list(list_on);
00243 cube_off = list_make_cube(list_off, noff);
00244 if (cube_off == NULL) {
00245 cpl_msg_error(_id,"could not create data cube!");
00246 return -1;
00247 }
00248
00249 destroy_list(list_off);
00250 im_twi = average_with_rejection( cube_twi, cfg->loReject, cfg->hiReject );
00251 if (im_twi == NULL) {
00252 cpl_msg_error(_id, "average_with_rejection failed\n" );
00253 return -1;
00254 }
00255
00256 im_dark = average_with_rejection( cube_dark,
00257 cfg->loReject,
00258 cfg->hiReject );
00259
00260 if (im_dark == NULL) {
00261 cpl_msg_error (_id, "average_with_rejection failed\n" );
00262 return -1;
00263 }
00264
00265 im_on = average_with_rejection( cube_on, cfg->loReject, cfg->hiReject );
00266 if (im_on == NULL) {
00267 cpl_msg_error (_id," average_with_rejection failed\n" );
00268 return -1;
00269 }
00270
00271 im_off = average_with_rejection( cube_off, cfg->loReject, cfg->hiReject );
00272 if (im_off == NULL) {
00273 cpl_msg_error ( _id," average_with_rejection failed\n" );
00274 return -1;
00275 }
00276
00277 destroy_cube(cube_twi);
00278 destroy_cube(cube_on);
00279 destroy_cube(cube_dark);
00280 destroy_cube(cube_off);
00281
00282
00283
00284
00285 im_twi_sub = sub_image(im_twi, im_dark);
00286 if (im_twi_sub == NULL) {
00287 cpl_msg_error ( _id," sub_image failed\n" );
00288 return -1;
00289 }
00290
00291 im_on_sub = sub_image(im_on, im_off);
00292 if (im_on_sub == NULL) {
00293 cpl_msg_error ( _id," sub_image failed\n" );
00294 return -1;
00295 }
00296
00297
00298 destroy_image(im_dark);
00299 destroy_image(im_off);
00300 destroy_image(im_twi);
00301 destroy_image(im_on);
00302
00303
00304
00305
00306
00307
00308
00309
00310 mask_im = load_image(cfg->mask);
00311 if (mask_im == NULL) {
00312 cpl_msg_error ( _id," could not load static bad pixel mask\n" );
00313 return -1;
00314 }
00315
00316 im_twi_bad = multImageByMask(im_twi_sub, mask_im);
00317 if (im_twi_bad == NULL) {
00318 cpl_msg_error ( _id," could not carry out multImageByMask\n" );
00319 return -1;
00320 }
00321
00322 im_on_bad = multImageByMask(im_on_sub, mask_im);
00323 if (im_on_bad == NULL) {
00324 cpl_msg_error ( _id," could not carry out multImageByMask\n" );
00325 return -1;
00326 }
00327
00328 destroy_image(mask_im);
00329 destroy_image(im_on_sub);
00330 destroy_image(im_twi_sub);
00331
00332
00333
00334
00335
00336
00337
00338
00339 if (cfg->warpfixInd == 1) {
00340
00341
00342 im_twi_tilt = image_warp(im_twi_bad, cfg->kernel, cfg->polyFile);
00343 if (im_twi_tilt == NULL) {
00344 cpl_msg_error ( _id," could not carry out image_warp\n" );
00345 return -1;
00346 }
00347
00348 im_on_tilt = image_warp(im_on_bad, cfg->kernel,cfg->polyFile );
00349 if (im_on_tilt == NULL) {
00350 cpl_msg_error ( _id," could not carry out image_warp\n" );
00351 return -1;
00352 }
00353
00354 destroy_image(im_twi_bad);
00355 destroy_image(im_on_bad);
00356 }
00357
00358
00359
00360
00361
00362
00363
00364
00365 if (cfg->warpfixInd == 1) {
00366 im_twi = im_twi_tilt;
00367 im_on = im_on_tilt;
00368 } else {
00369 im_twi = im_twi_bad;
00370 im_on = im_on_bad;
00371 }
00372
00373 row_twi = cleanMeanOfColumns(im_twi, cfg->loReject, cfg->hiReject);
00374 if (row_twi == NULL) {
00375 cpl_msg_error ( _id," could not carry out cleanMeanOfColumns\n" );
00376 return -1;
00377 }
00378 row_on = cleanMeanOfColumns(im_on, cfg->loReject, cfg->hiReject);
00379
00380 if (row_on == NULL) {
00381 cpl_msg_error ( _id," could not carry out cleanMeanOfColumns\n" );
00382 return -1;
00383 }
00384
00385 row_on_vec = imageToVector(row_on);
00386 row_twi_vec = imageToVector(row_twi);
00387
00388
00389
00390
00391
00392
00393 im_on_norm = divImageByRow( im_on, row_on_vec );
00394 if (im_on_norm == NULL) {
00395 cpl_msg_error ( _id," could not carry out divImageByRow\n" );
00396 return -1;
00397 }
00398
00399
00400
00401
00402 flat_field = multRowToImage( im_on_norm, row_twi_vec );
00403
00404 if (flat_field == NULL) {
00405 cpl_msg_error ( _id," could not carry out divImageByRow\n" );
00406 return -1;
00407 }
00408
00409
00410 ff_img=cpl_image_wrap_float(flat_field->lx,flat_field->ly,flat_field->data);
00411
00412
00413
00414
00415 if(-1 == sinfoni_pro_dump_ima(ff_img,raw,sof,cfg->outName,
00416 PRO_MASTER_TWIFLAT,_id,config)) {
00417 cpl_msg_error(_id,"cannot dump ima %s", cfg->outName);
00418 for (i=0; i< cfg->nframes; i++) {
00419
00420 destroy_image(im[i]);
00421 }
00422
00423 cplist = cpl_propertylist_new() ;
00424 if ((cpl_error_code)((cplist = cpl_propertylist_load(cfg->outName, 0)) == NULL)) {
00425 cpl_msg_error(_id, "getting header from frame %s",cfg->outName);
00426 cpl_propertylist_delete(cplist) ;
00427 return -1 ;
00428 }
00429
00430 change_plist_twiflat(cplist, cfg->outName, cfg->nframes);
00431
00432 if (cpl_image_save(ff_img, cfg->outName, CPL_BPP_DEFAULT,
00433 cplist,CPL_IO_DEFAULT)!=CPL_ERROR_NONE) {
00434 cpl_msg_error(_id, "Cannot save the product %s",cfg->outName);
00435 cpl_propertylist_delete(cplist) ;
00436 return -1 ;
00437 }
00438 cpl_propertylist_delete(cplist) ;
00439
00440
00441
00442
00443
00444
00445
00446 destroy_image(flat_field);
00447 destroy_image(im_on_norm);
00448 destroy_image(im_on);
00449 destroy_image(im_twi);
00450
00451 destroyVector(row_twi_vec);
00452 destroyVector(row_on_vec);
00453 twiflat_free (cfg);
00454 return -1;
00455 }
00456
00457
00458 for (i=0; i< cfg->nframes; i++) {
00459
00460 destroy_image(im[i]);
00461 }
00462
00463
00464
00465 destroy_image(flat_field);
00466 destroy_image(im_on_norm);
00467 destroy_image(im_on);
00468 destroy_image(im_twi);
00469
00470 destroyVector(row_twi_vec);
00471 destroyVector(row_on_vec);
00472 twiflat_free (cfg);
00473
00474
00475
00476
00477
00478 return 0 ;
00479 }
00480
00481
00482 void change_plist_twiflat (cpl_propertylist* plist, char* out_name, int n_frames) {
00483 char* date;
00484 char* text;
00485
00486 char file_name[FILENAMESZ] ;
00487 char file_tmp[FILENAMESZ] ;
00488 text = cpl_calloc(MAX_NAME_SIZE,sizeof(char));
00489 date = get_datetime_iso8601() ;
00490
00491
00492 strcpy(file_tmp,out_name);
00493 strcpy(file_name,"twiflat");
00494 strcat(file_name," -f");
00495
00496
00497 if (strlen(file_name) > 30) {
00498 strcpy(file_name,"twiflat");
00499 strcat(file_name," -f ");
00500
00501 }
00502
00503
00504
00505 cpl_propertylist_insert_after_string(plist, "EXPTIME", KEY_NAME_PIPEFILE, out_name) ;
00506 cpl_propertylist_set_comment(plist, KEY_NAME_PIPEFILE,"pipeline result filename" ) ;
00507
00508 cpl_propertylist_insert_after_string(plist, "EXPTIME", KEY_NAME_HPRO_DID,
00509 "ESO-VLT-DIC.PRO-1.13") ;
00510 cpl_propertylist_set_comment(plist, KEY_NAME_HPRO_DID,"data dictionary version" ) ;
00511
00512 cpl_propertylist_insert_after_string(plist, "EXPTIME", KEY_NAME_HPRO_TYPE,
00513 "REDUCED") ;
00514 cpl_propertylist_set_comment(plist, KEY_NAME_HPRO_TYPE, "product type" ) ;
00515
00516 cpl_propertylist_insert_after_string(plist, "EXPTIME", KEY_NAME_HPRO_CATG,
00517 "MASTER_LAMP_TWIFLAT") ;
00518 cpl_propertylist_set_comment(plist, KEY_NAME_HPRO_CATG, "product category" ) ;
00519
00520 cpl_propertylist_insert_after_string(plist, "EXPTIME", KEY_NAME_HPRO_STATUS,
00521 "OK") ;
00522 cpl_propertylist_set_comment(plist, KEY_NAME_HPRO_STATUS, "pipeline status" ) ;
00523
00524 cpl_propertylist_insert_after_string(plist, "EXPTIME", KEY_NAME_HPRO_DATE,
00525 date) ;
00526 cpl_propertylist_set_comment(plist, KEY_NAME_HPRO_CATG, "pipeline execution date" ) ;
00527
00528 cpl_propertylist_insert_after_string(plist, "EXPTIME", KEY_NAME_HPRO_RECID,
00529 file_name) ;
00530 cpl_propertylist_set_comment(plist, KEY_NAME_HPRO_RECID, "recipe ID" ) ;
00531
00532 cpl_propertylist_insert_after_string(plist, "EXPTIME", KEY_NAME_HPRO_DRSID,
00533 PACKAGE_VERSION) ;
00534 cpl_propertylist_set_comment(plist, KEY_NAME_HPRO_DRSID, "data reduction system ID" ) ;
00535
00536 sprintf(text,"%d",n_frames);
00537 cpl_propertylist_insert_after_string(plist, "EXPTIME", KEY_NAME_HPRO_DATANCOM,
00538 text) ;
00539 cpl_propertylist_set_comment(plist, KEY_NAME_HPRO_DATANCOM, "# of combined frames" ) ;
00540
00541 cpl_free(text);
00542 }
00543