00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019 #ifdef HAVE_CONFIG_H
00020 # include <config.h>
00021 #endif
00022 #include "sinfo_utilities_scired.h"
00023 #include "sinfo_functions.h"
00024 #include "sinfo_pfits.h"
00025 #include "sinfo_utils_wrappers.h"
00033 int
00034 sinfo_check_input_data(object_config* cfg)
00035 {
00036
00037
00038 if (cfg == NULL)
00039 {
00040 sinfo_msg_error (" could not parse cpl input!\n") ;
00041 return -1 ;
00042 }
00043
00044
00045 if(sinfo_is_fits_file(cfg->wavemap) != 1) {
00046 sinfo_msg_error("Input file wavemap %s is not FITS",cfg->wavemap);
00047 return -1;
00048 }
00049
00050
00051 if (cfg->halocorrectInd == 1)
00052 {
00053 if(sinfo_is_fits_file(cfg->halospectrum) != 1) {
00054 sinfo_msg_error("Input file %s is not FITS",cfg->halospectrum);
00055 return -1;
00056 }
00057
00058 }
00059
00060 if (cfg->northsouthInd == 0) {
00061 if (sinfo_is_fits_file(cfg->poslist) != 1)
00062 {
00063 sinfo_msg_error("File %s with tag %s is not FITS!",
00064 cfg->poslist,PRO_SLIT_POS);
00065 return -1 ;
00066 }
00067 } else {
00068
00069 if (sinfo_is_fits_file(cfg->distlist) != 1)
00070 {
00071 sinfo_msg_error("File %s with tag %s is not FITS!",
00072 cfg->distlist,PRO_SLITLETS_DISTANCE);
00073 return -1;
00074 }
00075 }
00076
00077
00078 return 0;
00079
00080
00081 }
00082
00083
00084
00096 int
00097 sinfo_auto_size_cube3(char** framelist,
00098 const int nframes,
00099 float* ref_offx,
00100 float* ref_offy,
00101 int* size_x,
00102 int* size_y)
00103 {
00104
00105 char* name =NULL;
00106 int n=0;
00107 float offx=0;
00108 float offy=0;
00109 float min_offx=0;
00110 float max_offx=0;
00111 float min_offy=0;
00112 float max_offy=0;
00113 cpl_propertylist* plist=NULL;
00114 sinfo_msg ("Automatic computation of output cube size") ;
00115 for ( n = 0 ; n < nframes ; n++ ) {
00116 name = framelist[n] ;
00117 plist=cpl_propertylist_load(name,0);
00118 offx = sinfo_pfits_get_cumoffsetx(plist);
00119 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00120 sinfo_msg_warning(" could not read fits header keyword cummoffsetx!");
00121 sinfo_msg_warning(" set it to 0");
00122 offx = 0;
00123 cpl_error_reset();
00124
00125 }
00126
00127 offy = sinfo_pfits_get_cumoffsety(plist);
00128 sinfo_free_propertylist(&plist);
00129 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00130 sinfo_msg_warning(" could not read fits header keyword! cumoffsety");
00131 sinfo_msg_warning(" set it to 0");
00132 offy = 0;
00133 cpl_error_reset();
00134
00135 }
00136 if(n==0) {
00137 min_offx=offx;
00138 min_offy=offy;
00139 max_offx=offx;
00140 max_offy=offy;
00141 } else {
00142 if(offx > max_offx) max_offx=offx;
00143 if(offy > max_offy) max_offy=offy;
00144 if(offx < min_offx) min_offx=offx;
00145 if(offy < min_offy) min_offy=offy;
00146 }
00147 }
00148
00149 *ref_offx=(min_offx+max_offx)/2;
00150 *ref_offy=(min_offy+max_offy)/2;
00151
00152 if(*size_x == 0) *size_x=2*floor(max_offx-min_offx+0.5)+64;
00153 if(*size_y == 0) *size_y=2*floor(max_offy-min_offy+0.5)+64;
00154 sinfo_msg("Output cube size: %d x %d",*size_x,*size_y);
00155 sinfo_msg_debug("Ref offset. x: %f y: %f",*ref_offx,*ref_offy);
00156 sinfo_msg_debug("Max offset. x: %f y: %f",max_offx,max_offy);
00157 sinfo_msg_debug("Min offset. x: %f y: %f",min_offx,min_offy);
00158 return 0;
00159
00160 }
00161
00162
00163
00172 int
00173 sinfo_auto_size_cube2(object_config * cfg, float* ref_offx, float* ref_offy)
00174 {
00175
00176 char* name =NULL;
00177 int n=0;
00178 float offx=0;
00179 float offy=0;
00180 float min_offx=0;
00181 float max_offx=0;
00182 float min_offy=0;
00183 float max_offy=0;
00184 cpl_propertylist * plist=NULL;
00185 sinfo_msg ("Automatic computation of output cube size") ;
00186 for ( n = 0 ; n < cfg->nframes ; n++ ) {
00187 name = cfg->framelist[n] ;
00188 plist=cpl_propertylist_load(name,0);
00189 offx = sinfo_pfits_get_cumoffsetx(plist);
00190 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00191 sinfo_msg_warning(" could not read fits header keyword cummoffsetx!");
00192 sinfo_msg_warning(" set it to 0");
00193 offx = 0;
00194 cpl_error_reset();
00195
00196 }
00197
00198 offy = sinfo_pfits_get_cumoffsety(plist);
00199 sinfo_free_propertylist(&plist);
00200 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00201 sinfo_msg_warning(" could not read fits header keyword! cumoffsety") ;
00202 sinfo_msg_warning(" set it to 0");
00203 offy = 0;
00204 cpl_error_reset();
00205
00206 }
00207 if(n==0) {
00208 min_offx=offx;
00209 min_offy=offy;
00210 max_offx=offx;
00211 max_offy=offy;
00212 } else {
00213 if(offx > max_offx) max_offx=offx;
00214 if(offy > max_offy) max_offy=offy;
00215 if(offx < min_offx) min_offx=offx;
00216 if(offy < min_offy) min_offy=offy;
00217 }
00218 }
00219 *ref_offx=(min_offx+max_offx)/2;
00220 *ref_offy=(min_offy+max_offy)/2;
00221
00222 if(cfg->size_x == 0) cfg->size_x=2*floor(max_offx-min_offx+0.5)+64;
00223 if(cfg->size_y == 0) cfg->size_y=2*floor(max_offy-min_offy+0.5)+64;
00224 sinfo_msg("Output cube size: %d x %d",cfg->size_x,cfg->size_y);
00225 sinfo_msg("Ref offset. x: %f y: %f",*ref_offx,*ref_offy);
00226 sinfo_msg("Max offset. x: %f y: %f",max_offx,max_offy);
00227 sinfo_msg("Min offset. x: %f y: %f",min_offx,min_offy);
00228 return 0;
00229
00230
00231 }
00232
00241 float*
00242 sinfo_read_distances(object_config* cfg)
00243 {
00244 int i=0;
00245 int* status=NULL;
00246 float * distances = NULL;
00247 float tmp_float=0;
00248 char tbl_distances_name[FILE_NAME_SZ];
00249 cpl_table* tbl_distances = NULL;
00250
00251 sinfo_msg("Read distances");
00252 distances = (float*) cpl_calloc (cfg->nslits - 1, sizeof (float));
00253
00254 if ( NULL == distances )
00255 {
00256 sinfo_msg_error ("could allocate memory!") ;
00257 return NULL ;
00258 }
00259
00260
00261 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00262 sinfo_msg_error("Before loading input table");
00263 sinfo_msg_error((char* ) cpl_error_get_message());
00264 return NULL;
00265 }
00266 strcpy(tbl_distances_name,cfg->distlist);
00267 tbl_distances = cpl_table_load(tbl_distances_name,1,0);
00268 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00269 sinfo_msg_error("loading input table %s",tbl_distances_name);
00270 sinfo_msg_error((char* ) cpl_error_get_message());
00271 return NULL;
00272 }
00273
00274 for (i =0 ; i< cfg->nslits-1; i++){
00275 tmp_float=cpl_table_get_float(tbl_distances,"slitlet_distance",i,status);
00276 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00277 sinfo_msg_error("reading col %s from table %s","slitlet_distance",
00278 tbl_distances_name);
00279 sinfo_msg_error((char* ) cpl_error_get_message());
00280 return NULL;
00281 }
00282 sinfo_new_array_set_value(distances,tmp_float,i);
00283 }
00284 cpl_table_delete(tbl_distances);
00285 return distances;
00286
00287 }
00288
00289
00290
00297 float**
00298 sinfo_read_slitlets_edges(object_config * cfg)
00299 {
00300
00301
00302 char tbl_slitpos_name[FILE_NAME_SZ];
00303 cpl_table* tbl_slitpos=NULL;
00304 int n=0;
00305 int i=0;
00306 int* status=NULL;
00307 float edge_x=0;
00308 float edge_y=0;
00309 float ** slit_edges = NULL;
00310
00311 slit_edges = sinfo_new_2Dfloatarray(cfg->nslits, 2) ;
00312
00313 strcpy(tbl_slitpos_name,cfg->poslist);
00314 tbl_slitpos = cpl_table_load(tbl_slitpos_name,1,0);
00315 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00316 sinfo_msg_error("error loading tbl %s",tbl_slitpos_name);
00317 sinfo_msg_error((char* ) cpl_error_get_message());
00318 return NULL;
00319 }
00320 n = cpl_table_get_nrow(tbl_slitpos);
00321 if (n != cfg->nslits) {
00322 sinfo_msg_error("No of slitlets in table is n = %d != %d !",n,cfg->nslits);
00323 return NULL;
00324 }
00325
00326 for (i =0 ; i< cfg->nslits; i++){
00327 edge_x=cpl_table_get_double(tbl_slitpos,"pos1",i,status);
00328 edge_y=cpl_table_get_double(tbl_slitpos,"pos2",i,status);
00329 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00330 sinfo_msg_error("error reading tbl %s row %d",tbl_slitpos_name,i);
00331 sinfo_msg_error((char* ) cpl_error_get_message());
00332 return NULL;
00333 }
00334 sinfo_new_array2D_set_value(slit_edges,edge_x,i,0);
00335 sinfo_new_array2D_set_value(slit_edges,edge_y,i,1);
00336 }
00337 cpl_table_delete(tbl_slitpos);
00338 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00339 sinfo_msg_error("error reading tbl %s",tbl_slitpos_name);
00340 sinfo_msg_error((char* ) cpl_error_get_message());
00341 return NULL;
00342 }
00343
00344 return slit_edges;
00345
00346 }
00347
00348
00349
00372
00373 cpl_imagelist * sinfo_new_cube_getvig(
00374 cpl_imagelist * cube_in,
00375 int loleft_x,
00376 int loleft_y,
00377 int upright_x,
00378 int upright_y)
00379 {
00380 cpl_imagelist * cube_out ;
00381 int i ;
00382 int outlx,
00383 outly ;
00384
00385 int ilx=0;
00386 int ily=0;
00387 int inp=0;
00388 cpl_image* i_img=NULL;
00389 cpl_image* o_img=NULL;
00390
00391 if (cube_in==NULL) return NULL ;
00392
00393 ilx=cpl_image_get_size_x(cpl_imagelist_get(cube_in,0));
00394 ily=cpl_image_get_size_y(cpl_imagelist_get(cube_in,0));
00395 inp=cpl_imagelist_get_size(cube_in);
00396
00397 if ((loleft_x>upright_x) ||
00398 (loleft_y>upright_y)) {
00399 sinfo_msg_error("ill-defined slit for extraction: aborting");
00400 return NULL ;
00401 }
00402
00403
00404 outlx = upright_x - loleft_x + 1 ;
00405 outly = upright_y - loleft_y + 1 ;
00406
00407
00408
00409 cube_out = cpl_imagelist_new() ;
00410
00411 for (i=0 ; i<cpl_imagelist_get_size(cube_in) ; i++) {
00412
00413 i_img=cpl_imagelist_get(cube_in,i);
00414
00415 o_img = sinfo_new_image_getvig(i_img,
00416 loleft_x, loleft_y,
00417 upright_x, upright_y) ;
00418 cpl_imagelist_set(cube_out,o_img,i);
00419 }
00420 return cube_out ;
00421 }
00422
00439
00440
00441 cpl_image * sinfo_new_image_getvig(
00442 cpl_image * image_in,
00443 int loleft_x,
00444 int loleft_y,
00445 int upright_x,
00446 int upright_y)
00447 {
00448 cpl_image * slit_img ;
00449 int i, j ;
00450 register
00451 pixelvalue * inpt,
00452 * outpt ;
00453 int outlx, outly ;
00454 int ilx=0;
00455 int ily=0;
00456
00457 float* pidata=NULL;
00458 float* podata=NULL;
00459
00460 if (image_in==NULL) return NULL ;
00461
00462 ilx=cpl_image_get_size_x(image_in);
00463 ily=cpl_image_get_size_y(image_in);
00464 pidata=cpl_image_get_data_float(image_in);
00465
00466 if ((loleft_x<1) || (loleft_x>ilx) ||
00467 (loleft_y<1) || (loleft_y>ily) ||
00468 (upright_x<1) || (upright_x>ilx) ||
00469 (upright_y<1) || (upright_y>ily) ||
00470 (loleft_x>upright_x) || (loleft_y>upright_y)) {
00471 sinfo_msg_error("extraction zone is [%d %d] [%d %d]\n"
00472 "cannot extract such zone: aborting slit extraction",
00473 loleft_x, loleft_y, upright_x, upright_y) ;
00474 return NULL ;
00475 }
00476
00477 outlx = upright_x - loleft_x + 1 ;
00478 outly = upright_y - loleft_y + 1 ;
00479 slit_img = cpl_image_new(outlx, outly,CPL_TYPE_FLOAT) ;
00480 podata=cpl_image_get_data_float(slit_img);
00481
00482 for (j=0 ; j<outly ; j++) {
00483 inpt = pidata+loleft_x-1 + (j+loleft_y-1)*ilx ;
00484 outpt = podata + j*outlx ;
00485 for (i=0 ; i<outlx ; i++) {
00486 *outpt++ = *inpt++ ;
00487 }
00488 }
00489 return slit_img ;
00490 }
00491
00503 int
00504 sinfo_new_set_wcs_cube(cpl_imagelist* cub, const char* name, double clambda,
00505 double dis, double cpix, double cx, double cy)
00506 {
00507 cpl_propertylist* plist=NULL;
00508
00509 if ((cpl_error_code)((plist = cpl_propertylist_load(name, 0)) == NULL)) {
00510 sinfo_msg_error( "getting header from frame %s",name);
00511 cpl_propertylist_delete(plist) ;
00512 return -1 ;
00513 }
00514
00515 sinfo_new_change_plist_cube(plist, clambda, dis, cpix, cx, cy) ;
00516
00517 if (cpl_imagelist_save(cub, name, CPL_BPP_DEFAULT,
00518 plist,CPL_IO_DEFAULT)!=CPL_ERROR_NONE) {
00519 sinfo_msg_error( "Cannot save the product %s",name);
00520 cpl_propertylist_delete(plist) ;
00521 return -1 ;
00522 }
00523 cpl_propertylist_delete(plist) ;
00524 return 0;
00525
00526 }
00535 int
00536 sinfo_new_set_wcs_image(cpl_image* img,
00537 const char* name,
00538 double cx,
00539 double cy)
00540 {
00541 cpl_propertylist* plist=NULL;
00542
00543 if ((cpl_error_code)((plist = cpl_propertylist_load(name, 0)) == NULL)) {
00544 sinfo_msg_error( "getting header from frame %s",name);
00545 cpl_propertylist_delete(plist) ;
00546 return -1 ;
00547 }
00548 sinfo_new_change_plist_image(plist, cx, cy) ;
00549
00550 if (cpl_image_save(img, name, CPL_BPP_DEFAULT,
00551 plist,CPL_IO_DEFAULT)!=CPL_ERROR_NONE) {
00552 sinfo_msg_error( "Cannot save the product %s",name);
00553 cpl_propertylist_delete(plist) ;
00554 return -1 ;
00555 }
00556 cpl_propertylist_delete(plist) ;
00557 return 0;
00558
00559 }
00560
00561
00571 int
00572 sinfo_new_set_wcs_spectrum(cpl_image* img, const char* name, double clambda,
00573 double dis, double cpix)
00574 {
00575 cpl_propertylist* plist=NULL;
00576
00577 if ((cpl_error_code)((plist = cpl_propertylist_load(name, 0)) == NULL)) {
00578 sinfo_msg_error( "getting header from frame %s",name);
00579 cpl_propertylist_delete(plist) ;
00580 return -1 ;
00581 }
00582 sinfo_new_change_plist_spectrum(plist, clambda, dis,cpix) ;
00583
00584 if (cpl_image_save(img, name, CPL_BPP_DEFAULT,
00585 plist,CPL_IO_DEFAULT)!=CPL_ERROR_NONE) {
00586 sinfo_msg_error( "Cannot save the product %s",name);
00587 cpl_propertylist_delete(plist) ;
00588 return -1 ;
00589 }
00590 cpl_propertylist_delete(plist) ;
00591 return 0;
00592
00593 }
00605 void sinfo_new_change_plist_cube (cpl_propertylist * plist,
00606 float cenLambda,
00607 float dispersion,
00608 int cenpix,
00609 float center_x,
00610 float center_y )
00611 {
00612
00613 float pixelscale ;
00614 double ra ;
00615 double dec ;
00616 double angle ;
00617 float radangle ;
00618 float cd1_1, cd1_2, cd2_1, cd2_2 ;
00619 char firsttext[2*FILE_NAME_SZ] ;
00620 strcpy(firsttext, "sinfo_rec_objnod -f \0") ;
00621
00622 pixelscale = sinfo_pfits_get_pixscale(plist)/2. ;
00623 ra = sinfo_pfits_get_ra(plist) ;
00624 dec = sinfo_pfits_get_DEC(plist) ;
00625 angle = sinfo_pfits_get_posangle(plist) ;
00626
00627 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00628 cpl_error_reset();
00629 }
00630
00631
00632
00633 radangle = angle * PI_NUMB / 180. ;
00634 cd1_1 = cos(radangle) * pixelscale / 3600. ;
00635 cd1_2 = sin(radangle) * pixelscale / 3600. ;
00636 cd2_1 = -sin(radangle) * pixelscale / 3600. ;
00637 cd2_2 = cos(radangle) * pixelscale / 3600. ;
00638
00639 cpl_propertylist_erase_regexp(plist, "^CTYPE1",0);
00640 cpl_propertylist_insert_after_string(plist,"EXPTIME","CTYPE1","RA---TAN");
00641 cpl_propertylist_set_comment(plist, "CTYPE1", "Projected Rectascension");
00642 cpl_propertylist_erase_regexp(plist, "^CRPIX1",0);
00643 cpl_propertylist_insert_after_float(plist,"CTYPE1","CRPIX1", center_x ) ;
00644 cpl_propertylist_set_comment(plist, "CRPIX1","Reference pixel in RA" ) ;
00645
00646 cpl_propertylist_erase_regexp(plist, "^CRVAL1",0);
00647 cpl_propertylist_insert_after_double(plist, "CRPIX1", "CRVAL1", ra ) ;
00648 cpl_propertylist_set_comment(plist, "CRVAL1","Reference RA" ) ;
00649
00650 cpl_propertylist_erase_regexp(plist, "^CDELT1",0);
00651 cpl_propertylist_insert_after_float(plist,"CRVAL1","CDELT1",
00652 -pixelscale/3600. ) ;
00653 cpl_propertylist_set_comment(plist, "CDELT1","pixel scale" ) ;
00654
00655 cpl_propertylist_erase_regexp(plist, "^CUNIT1",0);
00656 cpl_propertylist_insert_after_string(plist, "CDELT1", "CUNIT1", "deg" ) ;
00657 cpl_propertylist_set_comment(plist, "CUNIT1","RA-UNIT" ) ;
00658
00659 cpl_propertylist_erase_regexp(plist, "^CTYPE2",0);
00660 cpl_propertylist_insert_after_string(plist,"CUNIT1","CTYPE2","DEC--TAN");
00661 cpl_propertylist_set_comment(plist, "CTYPE2", "Projected Declination") ;
00662
00663 cpl_propertylist_erase_regexp(plist, "^CRPIX2",0);
00664 cpl_propertylist_insert_after_float(plist, "CTYPE2", "CRPIX2", center_y ) ;
00665 cpl_propertylist_set_comment(plist, "CRPIX2", "Reference pixel in DEC") ;
00666
00667 cpl_propertylist_erase_regexp(plist, "^CRVAL2",0);
00668 cpl_propertylist_insert_after_double(plist, "CRPIX2", "CRVAL2", dec) ;
00669 cpl_propertylist_set_comment(plist, "CRVAL2", "Reference DEC") ;
00670
00671 cpl_propertylist_erase_regexp(plist, "^CDELT2",0);
00672 cpl_propertylist_insert_after_float(plist, "CRVAL2", "CDELT2",
00673 pixelscale/3600. ) ;
00674 cpl_propertylist_set_comment(plist, "CDELT2", "pixel scale") ;
00675
00676 cpl_propertylist_erase_regexp(plist, "^CUNIT2",0);
00677 cpl_propertylist_insert_after_string(plist, "CDELT2", "CUNIT2", "deg" ) ;
00678 cpl_propertylist_set_comment(plist, "CUNIT2", "DEC-UNIT") ;
00679
00680 cpl_propertylist_erase_regexp(plist, "^CTYPE3",0);
00681 cpl_propertylist_insert_after_string(plist, "CUNIT2", "CTYPE3", "WAVE" ) ;
00682 cpl_propertylist_set_comment(plist,"CTYPE3","wavelength axis in microns") ;
00683
00684 cpl_propertylist_erase_regexp(plist, "^CRPIX3",0);
00685 cpl_propertylist_insert_after_int(plist, "CTYPE3", "CRPIX3", cenpix ) ;
00686 cpl_propertylist_set_comment(plist, "CRPIX3", "Reference pixel in z") ;
00687
00688 cpl_propertylist_erase_regexp(plist, "^CRVAL3",0);
00689 cpl_propertylist_insert_after_float(plist,"CRPIX3", "CRVAL3", cenLambda) ;
00690 cpl_propertylist_set_comment(plist, "CRVAL3", "central wavelength") ;
00691
00692 cpl_propertylist_erase_regexp(plist, "^CDELT3",0);
00693
00694 cpl_propertylist_insert_after_float(plist,"CRVAL3","CDELT3",dispersion) ;
00695 cpl_propertylist_set_comment(plist, "CDELT3", "microns per pixel") ;
00696
00697 cpl_propertylist_erase_regexp(plist, "^CUNIT3",0);
00698 cpl_propertylist_insert_after_string(plist,"CDELT3", "CUNIT3", "mum" ) ;
00699 cpl_propertylist_set_comment(plist, "CUNIT3", "spectral unit" ) ;
00700
00701 cpl_propertylist_erase_regexp(plist, "^CD1_1",0);
00702 cpl_propertylist_insert_after_float(plist, "CUNIT3", "CD1_1", cd1_1 ) ;
00703 cpl_propertylist_set_comment(plist, "CD1_1", "CD rotation matrix" ) ;
00704
00705 cpl_propertylist_erase_regexp(plist, "^CD1_2",0);
00706 cpl_propertylist_insert_after_float(plist, "CD1_1", "CD1_2", cd1_2 ) ;
00707 cpl_propertylist_set_comment(plist, "CD1_2", "CD rotation matrix" ) ;
00708
00709 cpl_propertylist_erase_regexp(plist, "^CD2_1",0);
00710 cpl_propertylist_insert_after_float(plist, "CD1_2", "CD2_1", cd2_1 ) ;
00711 cpl_propertylist_set_comment(plist, "CD2_1", "CD rotation matrix" ) ;
00712
00713 cpl_propertylist_erase_regexp(plist, "^CD2_2",0);
00714 cpl_propertylist_insert_after_float(plist, "CD2_1", "CD2_2", cd2_2 ) ;
00715 cpl_propertylist_set_comment(plist, "CD2_2", "CD rotation matrix" ) ;
00716
00717
00718 }
00719
00728 void sinfo_new_change_plist_spectrum (cpl_propertylist * plist,
00729 double cenLambda,
00730 double dispersion,
00731 int cenpix)
00732 {
00733
00734 cpl_propertylist_erase_regexp(plist, "^CTYPE3",0);
00735 cpl_propertylist_erase_regexp(plist, "^CRPIX3",0);
00736 cpl_propertylist_erase_regexp(plist, "^CRVAL3",0);
00737 cpl_propertylist_erase_regexp(plist, "^CDELT3",0);
00738 cpl_propertylist_erase_regexp(plist, "^CUNIT3",0);
00739
00740 cpl_propertylist_erase_regexp(plist, "^CTYPE2",0);
00741 cpl_propertylist_erase_regexp(plist, "^CRPIX2",0);
00742 cpl_propertylist_erase_regexp(plist, "^CRVAL2",0);
00743 cpl_propertylist_erase_regexp(plist, "^CDELT2",0);
00744 cpl_propertylist_erase_regexp(plist, "^CUNIT2",0);
00745
00746
00747 cpl_propertylist_erase_regexp(plist, "^CD1_1",0);
00748 cpl_propertylist_erase_regexp(plist, "^CD1_2",0);
00749 cpl_propertylist_erase_regexp(plist, "^CD2_1",0);
00750 cpl_propertylist_erase_regexp(plist, "^CD2_2",0);
00751
00752 cpl_propertylist_erase_regexp(plist, "^CTYPE1",0);
00753 cpl_propertylist_insert_after_string(plist, "EXPTIME", "CTYPE1", "PIXEL");
00754 cpl_propertylist_set_comment(plist, "CTYPE1", "Pixel coordinate system.");
00755
00756 cpl_propertylist_erase_regexp(plist, "^CRPIX1",0);
00757 cpl_propertylist_insert_after_int(plist, "CTYPE1", "CRPIX1", 1 ) ;
00758 cpl_propertylist_set_comment(plist, "CRPIX1", "Reference pixel in x") ;
00759
00760 cpl_propertylist_erase_regexp(plist, "^CRVAL1",0);
00761 cpl_propertylist_insert_after_double(plist, "CRPIX1", "CRVAL1", 1 ) ;
00762 cpl_propertylist_set_comment(plist, "CRVAL1", "value of ref pixel.") ;
00763
00764 cpl_propertylist_erase_regexp(plist, "^CDELT1",0);
00765 cpl_propertylist_insert_after_double(plist, "CRVAL1", "CDELT1", 1 ) ;
00766 cpl_propertylist_set_comment(plist, "CDELT1", "pixel scale") ;
00767
00768
00769
00770
00771
00772 cpl_propertylist_erase_regexp(plist, "^CUNIT1",0);
00773 cpl_propertylist_insert_after_string(plist, "CDELT1", "CUNIT1", "Pixel" );
00774 cpl_propertylist_set_comment(plist, "CUNIT1", "spectral unit" );
00775
00776 cpl_propertylist_erase_regexp(plist, "^CTYPE2",0);
00777 cpl_propertylist_insert_after_string(plist, "EXPTIME","CTYPE2","WAVE" );
00778 cpl_propertylist_set_comment(plist, "CTYPE2","wavelength axis in microns");
00779
00780 cpl_propertylist_erase_regexp(plist, "^CRPIX2",0);
00781 cpl_propertylist_insert_after_int(plist, "CTYPE2", "CRPIX2",cenpix ) ;
00782 cpl_propertylist_set_comment(plist, "CRPIX2", "Reference pixel in x") ;
00783
00784 cpl_propertylist_erase_regexp(plist, "^CRVAL2",0);
00785 cpl_propertylist_insert_after_double(plist, "CRPIX2","CRVAL2",cenLambda ) ;
00786 cpl_propertylist_set_comment(plist, "CRVAL2", "central wavelength") ;
00787
00788 cpl_propertylist_erase_regexp(plist, "^CDELT2",0);
00789 cpl_propertylist_insert_after_double(plist, "CRVAL2", "CDELT2",dispersion);
00790 cpl_propertylist_set_comment(plist, "CDELT2", "microns per pixel");
00791
00792 cpl_propertylist_erase_regexp(plist, "^CUNIT2",0);
00793 cpl_propertylist_insert_after_string(plist, "CDELT2", "CUNIT2", "mum");
00794 cpl_propertylist_set_comment(plist, "CUNIT2", "spectral unit" );
00795
00796
00797
00798 }
00808 void sinfo_new_change_plist_image (cpl_propertylist * plist,
00809 float center_x,
00810 float center_y )
00811 {
00812
00813 float pixelscale ;
00814 double ra ;
00815 double dec ;
00816 double angle ;
00817 float radangle ;
00818 float cd1_1, cd1_2, cd2_1, cd2_2 ;
00819 char firsttext[2*FILE_NAME_SZ] ;
00820
00821 strcpy(firsttext, "sinfo_rec_objnod -f \0") ;
00822
00823 pixelscale = sinfo_pfits_get_pixscale(plist)/2. ;
00824 ra = sinfo_pfits_get_ra(plist) ;
00825 dec = sinfo_pfits_get_DEC(plist) ;
00826 angle = sinfo_pfits_get_posangle(plist) ;
00827
00828 if(cpl_error_get_code() != CPL_ERROR_NONE) {
00829 cpl_error_reset();
00830 }
00831
00832 radangle = angle * PI_NUMB / 180. ;
00833 cd1_1 = cos(radangle) * pixelscale / 3600. ;
00834 cd1_2 = sin(radangle) * pixelscale / 3600. ;
00835 cd2_1 = -sin(radangle) * pixelscale / 3600. ;
00836 cd2_2 = cos(radangle) * pixelscale / 3600. ;
00837
00838 cpl_propertylist_erase_regexp(plist, "^CTYPE1",0);
00839 cpl_propertylist_insert_after_string(plist,"EXPTIME","CTYPE1","RA---TAN") ;
00840 cpl_propertylist_set_comment(plist, "CTYPE1", "Projected Rectascension" ) ;
00841
00842 cpl_propertylist_erase_regexp(plist, "^CRPIX1",0);
00843 cpl_propertylist_insert_after_float(plist, "CTYPE1", "CRPIX1", center_x) ;
00844 cpl_propertylist_set_comment(plist, "CRPIX1","Reference pixel in RA" ) ;
00845
00846 cpl_propertylist_erase_regexp(plist, "^CRVAL1",0);
00847 cpl_propertylist_insert_after_double(plist, "CRPIX1", "CRVAL1", ra ) ;
00848 cpl_propertylist_set_comment(plist, "CRVAL1","Reference RA" ) ;
00849
00850 cpl_propertylist_erase_regexp(plist, "^CDELT1",0);
00851 cpl_propertylist_insert_after_float(plist,"CRVAL1","CDELT1",
00852 -pixelscale/3600. ) ;
00853 cpl_propertylist_set_comment(plist, "CDELT1","pixel scale" ) ;
00854
00855 cpl_propertylist_erase_regexp(plist, "^CUNIT1",0);
00856 cpl_propertylist_insert_after_string(plist, "CDELT1", "CUNIT1", "deg" ) ;
00857 cpl_propertylist_set_comment(plist, "CUNIT1","RA-UNIT" ) ;
00858
00859 cpl_propertylist_erase_regexp(plist, "^CTYPE2",0);
00860 cpl_propertylist_insert_after_string(plist,"CUNIT1","CTYPE2","DEC--TAN") ;
00861 cpl_propertylist_set_comment(plist, "CTYPE2", "Projected Declination") ;
00862
00863 cpl_propertylist_erase_regexp(plist, "^CRPIX2",0);
00864 cpl_propertylist_insert_after_float(plist, "CTYPE2", "CRPIX2", center_y) ;
00865 cpl_propertylist_set_comment(plist, "CRPIX2", "Reference pixel in DEC") ;
00866
00867 cpl_propertylist_erase_regexp(plist, "^CRVAL2",0);
00868 cpl_propertylist_insert_after_double(plist, "CRPIX2", "CRVAL2", dec) ;
00869 cpl_propertylist_set_comment(plist, "CRVAL2", "Reference DEC") ;
00870
00871 cpl_propertylist_erase_regexp(plist, "^CDELT2",0);
00872 cpl_propertylist_insert_after_float(plist, "CRVAL2", "CDELT2",
00873 pixelscale/3600. ) ;
00874 cpl_propertylist_set_comment(plist, "CDELT2", "pixel scale") ;
00875
00876 cpl_propertylist_erase_regexp(plist, "^CUNIT2",0);
00877 cpl_propertylist_insert_after_string(plist, "CDELT2", "CUNIT2", "deg" ) ;
00878 cpl_propertylist_set_comment(plist, "CUNIT2", "DEC-UNIT") ;
00879
00880 cpl_propertylist_erase_regexp(plist, "^CD1_1",0);
00881 cpl_propertylist_insert_after_float(plist, "CUNIT2", "CD1_1", cd1_1 ) ;
00882 cpl_propertylist_set_comment(plist, "CD1_1", "CD rotation matrix" ) ;
00883
00884 cpl_propertylist_erase_regexp(plist, "^CD1_2",0);
00885 cpl_propertylist_insert_after_float(plist, "CD1_1", "CD1_2", cd1_2 ) ;
00886 cpl_propertylist_set_comment(plist, "CD1_2", "CD rotation matrix" ) ;
00887
00888 cpl_propertylist_erase_regexp(plist, "^CD2_1",0);
00889 cpl_propertylist_insert_after_float(plist, "CD1_2", "CD2_1", cd2_1 ) ;
00890 cpl_propertylist_set_comment(plist, "CD2_1", "CD rotation matrix" ) ;
00891
00892 cpl_propertylist_erase_regexp(plist, "^CD2_2",0);
00893 cpl_propertylist_insert_after_float(plist, "CD2_1", "CD2_2", cd2_2 ) ;
00894 cpl_propertylist_set_comment(plist, "CD2_2", "CD rotation matrix" ) ;
00895
00896
00897 }
00898
00906 cpl_imagelist**
00907 sinfo_new_sinfoni_correct_median(cpl_imagelist** cubes, const int n_cubes)
00908 {
00909 int i=0;
00910 cpl_imagelist** cubes_cor=NULL;
00911 double local_median=0;
00912 int z=0;
00913 cpl_image* i_img=NULL;
00914 cpl_image* o_img=NULL;
00915
00916
00917 if ( cubes == NULL ) {
00918 sinfo_msg_error ("no cube list given!") ;
00919 return NULL ;
00920 }
00921 if ( n_cubes <= 0 ) {
00922 sinfo_msg_error ("wrong number of data cubes in list!") ;
00923 return NULL ;
00924 }
00925
00926 cubes_cor = (cpl_imagelist**) cpl_calloc (n_cubes, sizeof (cpl_imagelist*)) ;
00927
00928 for ( i = 0 ; i < n_cubes ; i++ ) {
00929 cubes_cor[i] = cpl_imagelist_new();
00930 for(z=0;z< cpl_imagelist_get_size(cubes[i]); z++) {
00931 i_img=cpl_imagelist_get(cubes[i],z);
00932 local_median=cpl_image_get_median(i_img);;
00933 o_img=cpl_image_duplicate(i_img);
00934 if(!isnan(local_median)) {
00935 cpl_image_subtract_scalar(o_img,local_median);
00936 }
00937 cpl_imagelist_set(cubes_cor[i],o_img,z);
00938 }
00939 }
00940
00941 return cubes_cor;
00942 }
00943
00950 int sinfo_new_sinfoni_correct_median_it(cpl_imagelist** inp)
00951 {
00952
00953 double local_median=0;
00954 int z=0;
00955 cpl_image* img=NULL;
00956
00957 for(z=0;z< cpl_imagelist_get_size((*inp)); z++) {
00958 img=cpl_imagelist_get((*inp),z);
00959 local_median=sinfo_new_my_median_image(img);
00960 if(!isnan(local_median)) {
00961 cpl_image_subtract_scalar(img,local_median);
00962 } else {
00963 sinfo_msg_error("local_median is NAN");
00964 }
00965 cpl_imagelist_set((*inp),img,z);
00966 }
00967
00968 return 0;
00969 }
00970
00980 cpl_imagelist** sinfo_new_sinfoni_correct_sky(cpl_imagelist** cubes,
00981 const int nc,
00982 cpl_imagelist* sky_cube)
00983
00984 {
00985 cpl_imagelist** cubes_sky=NULL;
00986 int x=0;
00987 int y=0;
00988 int z=0;
00989 int i=0;
00990 float k=0.5;
00991 int ovr=0;
00992 int ks=0;
00993 int nclip=0;
00994 double med=0;
00995 double avg=0;
00996 double sig=0;
00997 int msk_sum=0;
00998 double val_msk_sum=0;
00999 cpl_vector* val=NULL;
01000 cpl_vector* msk=NULL;
01001 int ilx=0;
01002 int ily=0;
01003 int inp=0;
01004
01005
01006 int sky_lx=0;
01007 int sky_ly=0;
01008
01009 float* pidata=NULL;
01010 float* p_sky_data=NULL;
01011 cpl_image* i_img=NULL;
01012 cpl_image* sky_img=NULL;
01013
01014
01015 if ( cubes == NULL ) {
01016 sinfo_msg_error ("no cube list given!") ;
01017 return NULL ;
01018 }
01019 ilx=cpl_image_get_size_x(cpl_imagelist_get(cubes[0],0));
01020 ily=cpl_image_get_size_y(cpl_imagelist_get(cubes[0],0));
01021 inp=cpl_imagelist_get_size(cubes[0]);
01022
01023
01024 sky_lx=cpl_image_get_size_x(cpl_imagelist_get(sky_cube,0));
01025 sky_ly=cpl_image_get_size_y(cpl_imagelist_get(sky_cube,0));
01026 if ( nc <= 0 ) {
01027 sinfo_msg_error ("wrong number of data cubes in list!") ;
01028 return NULL ;
01029 }
01030
01031 cubes_sky = (cpl_imagelist**) cpl_calloc (nc, sizeof (cpl_imagelist*)) ;
01032 for(z=0;z< inp; z++) {
01033 sky_img=cpl_imagelist_get(sky_cube,z);
01034 p_sky_data=cpl_image_get_data_float(sky_img);
01035 for(y=0;y< ily; y++) {
01036 for(x=0;x< ilx; x++) {
01037
01038 msk=cpl_vector_new(nc);
01039 for (i=0;i<nc;i++) {
01040 cpl_vector_set(msk,i,1);
01041 }
01042 nclip=0;
01043 for (ks=0;ks<nc;ks++) {
01044 sig=0;
01045 med=0;
01046 ovr=0;
01047 val=cpl_vector_new(nc-nclip);
01048
01049 for ( i = 0 ; i < nc ; i++ ) {
01050 i_img=cpl_imagelist_get(cubes[i],z);
01051 pidata=cpl_image_get_data_float(i_img);
01052 if ((!isnan(pidata[x+y*ilx])) &&
01053 (cpl_vector_get(msk,i) != 0) ) {
01054 cpl_vector_set(val,ovr,(double)pidata[x+y*ilx]);
01055 ovr++;
01056 }
01057 }
01058
01059 if(ovr>0) {
01060 avg=cpl_vector_get_mean(val);
01061 med=cpl_vector_get_median(val);
01062 if(ovr>1) {
01063 sig=cpl_vector_get_stdev(val);
01064 } else {
01065 sig=0;
01066 }
01067 } else {
01068 avg=cpl_vector_get(val,0);
01069 med=avg;
01070 sig=0;
01071 }
01072
01073 cpl_vector_delete(val);
01074 for ( i = 0 ; i < nc ; i++ ) {
01075 i_img=cpl_imagelist_get(cubes[i],z);
01076 pidata=cpl_image_get_data_float(i_img);
01077
01078 if ((!isnan(pidata[x+y*ilx])) &&
01079 (cpl_vector_get(msk,i) != 0)) {
01080 if(abs((pidata[x+y*ilx]-med))> k*sig) {
01081 cpl_vector_set(msk,i,0);
01082 nclip++;
01083 }
01084 }
01085 }
01086 }
01087 msk_sum=0;
01088 val_msk_sum=0;
01089 for ( i = 0 ; i < nc ; i++ ) {
01090
01091 if (!isnan(pidata[x+y*ilx])) {
01092 msk_sum+=cpl_vector_get(msk,i);
01093 val_msk_sum+=pidata[x+y*ilx]*
01094 cpl_vector_get(msk,i);
01095 }
01096 }
01097 p_sky_data[x+y*sky_lx]=val_msk_sum/msk_sum;
01098 cpl_vector_delete(msk);
01099 }
01100 }
01101 }
01102 for ( i = 0 ; i < nc ; i++ ) {
01103 cubes_sky[i]=cpl_imagelist_duplicate(cubes[i]);
01104
01105 cpl_imagelist_subtract(cubes_sky[i],sky_cube);
01106
01107 }
01108
01109
01110 return cubes_sky;
01111 }
01112
01127 cpl_imagelist** sinfo_new_sinfoni_correct_sky2(cpl_imagelist** cubes,
01128 const int nc,
01129 cpl_imagelist* sky_cube,
01130 cpl_imagelist* med_cube,
01131 cpl_imagelist* msk_cube,
01132 cpl_imagelist* avg_cube,
01133 cpl_imagelist* sig_cube,
01134 cpl_imagelist* ovr_cube)
01135 {
01136 cpl_imagelist** cubes_sky=NULL;
01137 int x=0;
01138 int y=0;
01139 int z=0;
01140 int i=0;
01141 float k=0.5;
01142 int ovr=0;
01143 int ks=0;
01144 int nclip=0;
01145 double med=0;
01146 double avg=0;
01147 double sig=0;
01148 int msk_sum=0;
01149 double val_msk_sum=0;
01150 cpl_vector* val=NULL;
01151 cpl_vector* msk=NULL;
01152
01153 int ilx=0;
01154 int ily=0;
01155 int inp=0;
01156
01157 int ovr_lx=0;
01158 int msk_lx=0;
01159 int avg_lx=0;
01160 int sig_lx=0;
01161 int sky_lx=0;
01162 int med_lx=0;
01163
01164
01165 float* p_ovr_data=NULL;
01166 float* p_msk_data=NULL;
01167 float* p_avg_data=NULL;
01168 float* p_sig_data=NULL;
01169 float* p_sky_data=NULL;
01170 float* p_med_data=NULL;
01171
01172 float* pidata=NULL;
01173
01174 if ( cubes == NULL ) {
01175 sinfo_msg_error ("no cube list given!") ;
01176 return NULL ;
01177 }
01178
01179 ilx=cpl_image_get_size_x(cpl_imagelist_get(cubes[0],0));
01180 ily=cpl_image_get_size_y(cpl_imagelist_get(cubes[0],0));
01181 inp=cpl_imagelist_get_size(cubes[0]);
01182
01183 if ( nc <= 0 ) {
01184 sinfo_msg_error ("wrong number of data cubes in list!") ;
01185 return NULL ;
01186 }
01187
01188 cubes_sky = (cpl_imagelist**) cpl_calloc (nc, sizeof (cpl_imagelist*)) ;
01189
01190 ovr_lx=ilx;
01191 msk_lx=ilx;
01192 avg_lx=ilx;
01193 sig_lx=ilx;
01194 sky_lx=ilx;
01195 med_lx=ilx;
01196
01197 for(z=0;z< inp; z++) {
01198 p_ovr_data=cpl_image_get_data_float(cpl_imagelist_get(ovr_cube,z));
01199 p_msk_data=cpl_image_get_data_float(cpl_imagelist_get(msk_cube,z));
01200 p_avg_data=cpl_image_get_data_float(cpl_imagelist_get(avg_cube,z));
01201 p_sig_data=cpl_image_get_data_float(cpl_imagelist_get(sig_cube,z));
01202 p_sky_data=cpl_image_get_data_float(cpl_imagelist_get(sky_cube,z));
01203 p_med_data=cpl_image_get_data_float(cpl_imagelist_get(med_cube,z));
01204
01205
01206 for(y=0;y< ily; y++) {
01207 for(x=0;x< ilx; x++) {
01208
01209 msk=cpl_vector_new(nc);
01210 for (i=0;i<nc;i++) {
01211 cpl_vector_set(msk,i,1);
01212 }
01213 p_ovr_data[x+y*ovr_lx]=nc;
01214 p_msk_data[x+y*msk_lx]=nc;
01215 nclip=0;
01216 for (ks=0;ks<nc;ks++) {
01217 sig=0;
01218 med=0;
01219 ovr=0;
01220 val=cpl_vector_new(nc-nclip);
01221
01222 for ( i = 0 ; i < nc ; i++ ) {
01223 pidata=cpl_image_get_data_float(cpl_imagelist_get(cubes[i],z));
01224 if ((!isnan(pidata[x+y*ilx])) &&
01225 (cpl_vector_get(msk,i) != 0) ) {
01226 cpl_vector_set(val,ovr,(double)pidata[x+y*ilx]);
01227 ovr++;
01228 }
01229 }
01230
01231 if(ovr>1) {
01232 avg=cpl_vector_get_mean(val);
01233 med=cpl_vector_get_median(val);
01234 sig=cpl_vector_get_stdev(val);
01235 } else {
01236 avg=cpl_vector_get(val,0);
01237 med=avg;
01238 sig=0;
01239 }
01240
01241 p_med_data[x+y*med_lx]=med;
01242 p_avg_data[x+y*avg_lx]=avg;
01243 p_sig_data[x+y*sig_lx]=sig;
01244 cpl_vector_delete(val);
01245 for ( i = 0 ; i < nc ; i++ ) {
01246 pidata=cpl_image_get_data_float(cpl_imagelist_get(cubes[i],z));
01247
01248 if ((!isnan(pidata[x+y*ilx])) &&
01249 (cpl_vector_get(msk,i) != 0)) {
01250 if(abs((pidata[x+y*ilx]-med))> k*sig) {
01251
01252 p_msk_data[x+y*msk_lx]-=1;
01253 cpl_vector_set(msk,i,0);
01254 nclip++;
01255 }
01256 }
01257 }
01258 }
01259 msk_sum=0;
01260 val_msk_sum=0;
01261 for ( i = 0 ; i < nc ; i++ ) {
01262 pidata=cpl_image_get_data_float(cpl_imagelist_get(cubes[i],z));
01263
01264 if (!isnan(pidata[x+y*ilx])) {
01265
01266
01267
01268
01269
01270 msk_sum+=cpl_vector_get(msk,i);
01271 val_msk_sum+=pidata[x+y*ilx]*
01272 cpl_vector_get(msk,i);
01273 }
01274 }
01275 p_sky_data[x+y*sky_lx]=val_msk_sum/msk_sum;
01276 cpl_vector_delete(msk);
01277 }
01278 }
01279 }
01280 for ( i = 0 ; i < nc ; i++ ) {
01281 cubes_sky[i]=cpl_imagelist_duplicate(cubes[i]);
01282
01283 cpl_imagelist_subtract(cubes_sky[i],sky_cube);
01284 }
01285
01286
01287 return cubes_sky;
01288 }
01298 int
01299 sinfo_new_extract_ao_info(cpl_frameset* obj_set,
01300 cpl_frameset** set,
01301 const char* recid,
01302 cpl_parameterlist* config)
01303 {
01304
01305 cpl_frame* frm=NULL;
01306 cpl_table* tbl=NULL;
01307 cpl_frame* product_frame=NULL;
01308 cpl_propertylist* plist=NULL;
01309 cpl_image* image=NULL;
01310
01311 int nobj=0;
01312 int i=0;
01313 const char* name_o="out_ao_info.fits";
01314
01315 name_o="out_ao_info.fits";
01316 nobj=cpl_frameset_get_size(obj_set);
01317 for(i=0;i<nobj;i++) {
01318 frm=cpl_frameset_get_frame(obj_set,i);
01319 tbl=cpl_table_load(cpl_frame_get_filename(frm),1,1);
01320 plist=cpl_propertylist_load(cpl_frame_get_filename(frm),0);
01321 if(i==0) {
01322 cpl_table_save(tbl,NULL,plist,name_o,CPL_IO_DEFAULT);
01323 } else {
01324 cpl_table_save(tbl,NULL,plist,name_o,CPL_IO_EXTEND);
01325 }
01326 }
01327 image=cpl_image_new(0,0,CPL_TYPE_FLOAT);
01328
01329 product_frame = cpl_frame_new();
01330 cpl_frame_set_filename(product_frame, name_o) ;
01331 cpl_frame_set_tag(product_frame, PRO_AO_INFO) ;
01332 cpl_frame_set_type(product_frame, CPL_FRAME_TYPE_IMAGE) ;
01333 cpl_frame_set_group(product_frame, CPL_FRAME_GROUP_PRODUCT) ;
01334 cpl_frame_set_level(product_frame, CPL_FRAME_LEVEL_FINAL) ;
01335 if (cpl_error_get_code()) {
01336 sinfo_msg_error( "Error while initialising the product frame") ;
01337 cpl_propertylist_delete(plist) ;
01338 cpl_frame_delete(product_frame) ;
01339 cpl_image_delete(image) ;
01340 return -1 ;
01341 }
01342
01343
01344 if (cpl_dfs_setup_product_header(plist, product_frame, *set, config,
01345 recid, "SINFONI", "?Dictionary?") != CPL_ERROR_NONE) {
01346 sinfo_msg_error( "Problem in the product DFS-compliance") ;
01347 cpl_propertylist_delete(plist) ;
01348 cpl_frame_delete(product_frame) ;
01349 cpl_image_delete(image) ;
01350 return -1 ;
01351 }
01352
01353
01354 if (cpl_image_save(image, name_o, CPL_BPP_DEFAULT, plist,
01355 CPL_IO_DEFAULT) != CPL_ERROR_NONE) {
01356 sinfo_msg_error( "Could not save product");
01357 cpl_propertylist_delete(plist) ;
01358 cpl_frame_delete(product_frame) ;
01359 cpl_image_delete(image) ;
01360 return -1 ;
01361 }
01362 cpl_propertylist_delete(plist) ;
01363 cpl_image_delete(image) ;
01364
01365
01366 cpl_frameset_insert(*set, product_frame) ;
01367
01368 return 0;
01369
01370 }
01383 int
01384 sinfo_new_assign_offset(const int n,
01385 const char* name,
01386 float* offsetx,
01387 float* offsety,
01388 const float ref_offx,
01389 const float ref_offy)
01390 {
01391
01392 float offx=0;
01393 float offy=0;
01394 double mjd_obs=0;
01395
01396
01397
01398
01399
01400
01401
01402
01403
01404
01405
01406
01407 cpl_propertylist * plist=NULL;
01408 sinfo_msg_debug("Assign offsets");
01409
01410 if ((cpl_error_code)((plist = cpl_propertylist_load(name, 0)) == NULL)) {
01411 sinfo_msg_error( "getting header from reference frame %s",name);
01412 cpl_propertylist_delete(plist) ;
01413 return -1 ;
01414 }
01415
01416 offx = sinfo_pfits_get_cumoffsetx(plist) - ref_offx ;
01417 if(cpl_error_get_code() != CPL_ERROR_NONE) {
01418 sinfo_msg_warning(" could not read fits header keyword cummoffsetx!") ;
01419 sinfo_msg_warning(" Set relative offset to 0 - %f!",ref_offx) ;
01420 offx = - ref_offx;
01421 cpl_error_reset();
01422
01423 }
01424
01425 offy = sinfo_pfits_get_cumoffsety(plist) - ref_offy ;
01426 if(cpl_error_get_code() != CPL_ERROR_NONE) {
01427 sinfo_msg_warning(" could not read fits header keyword! cumoffsety") ;
01428 sinfo_msg_warning(" Set relative offset to 0 - %f!",ref_offx) ;
01429 offy = - ref_offy;
01430 cpl_error_reset();
01431
01432 }
01433 sinfo_msg_debug("offx=%f offy=%f",offx,offy);
01434
01435 if (cpl_propertylist_contains(plist, KEY_NAME_MJD_OBS)) {
01436 mjd_obs=cpl_propertylist_get_double(plist, KEY_NAME_MJD_OBS);
01437 } else {
01438 sinfo_msg_error("keyword %s does not exist",KEY_NAME_MJD_OBS);
01439 cpl_propertylist_delete(plist) ;
01440 return -1;
01441 }
01442 cpl_propertylist_delete(plist) ;
01443
01444 if (mjd_obs > 53825. ) {
01445
01446 sinfo_msg("New cumoffset setting convention");
01447 sinfo_new_array_set_value(offsetx,2*offx,n);
01448 sinfo_new_array_set_value(offsety,2*offy,n);
01449 } else if ((mjd_obs > 53421.58210082 ) && (mjd_obs <= 53825.)){
01450
01451
01452
01453
01454
01455 sinfo_new_array_set_value(offsetx,-2*offx,n);
01456 sinfo_new_array_set_value(offsety,2*offy,n);
01457 } else {
01458
01459
01460
01461
01462
01463 sinfo_new_array_set_value(offsetx,2*offx,n);
01464 sinfo_new_array_set_value(offsety,-2*offy,n);
01465 }
01466
01467 return 0;
01468
01469
01470 }
01483 int
01484 sinfo_new_object_assign_offset(const char* name,
01485 const int n,
01486 double* ref_offx,
01487 double* ref_offy,
01488 float** offsetx,
01489 float** offsety)
01490 {
01491
01492 float offx=0;
01493 float offy=0;
01494 double mjd_obs=0;
01495 cpl_propertylist * plist=NULL;
01496 sinfo_msg_debug("Assign offsets");
01497
01498 if ((cpl_error_code)((plist = cpl_propertylist_load(name, 0)) == NULL)) {
01499 sinfo_msg_error( "getting header from reference frame %s",name);
01500 cpl_propertylist_delete(plist) ;
01501 return -1 ;
01502 }
01503 if ( n == 0 ) {
01504
01505 *ref_offx = sinfo_pfits_get_cumoffsetx(plist) ;
01506 if(cpl_error_get_code() != CPL_ERROR_NONE) {
01507 sinfo_msg_error(" could not read fits header keyword cummoffsetx!") ;
01508
01509 cpl_error_reset();
01510 }
01511
01512 *ref_offy = sinfo_pfits_get_cumoffsety(plist) ;
01513 if(cpl_error_get_code() != CPL_ERROR_NONE) {
01514 sinfo_msg_error(" could not read fits header keyword! cumoffsety") ;
01515 cpl_error_reset();
01516
01517 }
01518 sinfo_msg_debug("Reference offx=%f offy=%f",*ref_offx,*ref_offy);
01519
01520 offx = 0. ;
01521 offy = 0. ;
01522
01523 } else {
01524
01525 offx = sinfo_pfits_get_cumoffsetx(plist) - *ref_offx ;
01526 if(cpl_error_get_code() != CPL_ERROR_NONE) {
01527 sinfo_msg_error(" could not read fits header keyword cummoffsetx!") ;
01528 cpl_error_reset();
01529
01530 }
01531
01532 offy = sinfo_pfits_get_cumoffsety(plist) - *ref_offy ;
01533 if(cpl_error_get_code() != CPL_ERROR_NONE) {
01534 sinfo_msg_error(" could not read fits header keyword! cumoffsety") ;
01535
01536 cpl_error_reset();
01537 }
01538 sinfo_msg_debug("offx=%f offy=%f",offx,offy);
01539 }
01540
01541
01542
01543
01544
01545
01546
01547
01548
01549
01550
01551
01552
01553
01554 if (cpl_propertylist_contains(plist, KEY_NAME_MJD_OBS)) {
01555 mjd_obs=cpl_propertylist_get_double(plist, KEY_NAME_MJD_OBS);
01556 } else {
01557 sinfo_msg_error("keyword %s does not exist",KEY_NAME_MJD_OBS);
01558 cpl_propertylist_delete(plist) ;
01559 return -1;
01560 }
01561 cpl_propertylist_delete(plist) ;
01562
01563 if (mjd_obs > 53825. ) {
01564
01565 sinfo_msg("New cumoffset setting convention");
01566 sinfo_new_array_set_value(*offsetx,2*offx,n);
01567 sinfo_new_array_set_value(*offsety,2*offy,n);
01568 } else if ((mjd_obs > 53421.58210082 ) && (mjd_obs <= 53825.)){
01569
01570 sinfo_new_array_set_value(*offsetx,-offx*2,n);
01571 sinfo_new_array_set_value(*offsety,+offy*2,n);
01572 } else {
01573
01574 sinfo_new_array_set_value(*offsetx,+offx*2,n);
01575 sinfo_new_array_set_value(*offsety,-offy*2,n);
01576 }
01577
01578 return 0;
01579 }
01580
01592 cpl_imagelist*
01593 sinfo_new_fine_tune(cpl_imagelist* cube,
01594 float* correct_dist,
01595 const char* method,
01596 const int order,
01597 const int nslits) {
01598 int i =0;
01599 cpl_imagelist* outcube2=NULL;
01600 float* neg_dist=NULL;
01601 sinfo_msg("Finetuning, method=%s",method);
01602
01603 if (strcmp(method,"P")==0)
01604 {
01605 outcube2 = sinfo_new_fine_tune_cube( cube, correct_dist, order ) ;
01606 if (outcube2 == NULL)
01607 {
01608 sinfo_msg_error (" could not fine tune the data cube\n") ;
01609 return NULL ;
01610 }
01611 }
01612 else if (strcmp(method,"F")==0)
01613 {
01614 neg_dist=cpl_calloc(nslits,sizeof(float));
01615 for ( i = 0 ; i < nslits ; i++ )
01616 {
01617 neg_dist[i] = -correct_dist[i] ;
01618 }
01619 outcube2 = sinfo_new_fine_tune_cube_by_FFT( cube, neg_dist ) ;
01620 cpl_free(neg_dist);
01621 if ( outcube2 == NULL )
01622 {
01623 sinfo_msg_error (" could not fine tune the data cube\n") ;
01624 return NULL ;
01625 }
01626 }
01627 else if (strcmp(method,"S")==0)
01628 {
01629 outcube2 = sinfo_new_fine_tune_cube_by_spline( cube, correct_dist ) ;
01630 if ( outcube2 == NULL )
01631 {
01632 sinfo_msg_error (" could not fine tune the data cube\n") ;
01633 return NULL ;
01634 }
01635 }
01636 else
01637 {
01638 sinfo_msg_error (" wrong method indicator given!") ;
01639 return NULL ;
01640 }
01641
01642
01643
01644 return outcube2;
01645
01646 }