KMOS Pipeline Reference Manual
1.3.11
|
00001 /* 00002 * This file is part of the KMOS Pipeline 00003 * Copyright (C) 2002,2003 European Southern Observatory 00004 * 00005 * This program is free software; you can redistribute it and/or modify 00006 * it under the terms of the GNU General Public License as published by 00007 * the Free Software Foundation; either version 2 of the License, or 00008 * (at your option) any later version. 00009 * 00010 * This program is distributed in the hope that it will be useful, 00011 * but WITHOUT ANY WARRANTY; without even the implied warranty of 00012 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00013 * GNU General Public License for more details. 00014 * 00015 * You should have received a copy of the GNU General Public License 00016 * along with this program; if not, write to the Free Software 00017 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 00018 */ 00019 00020 #ifdef HAVE_CONFIG_H 00021 #include <config.h> 00022 #endif 00023 00024 /*----------------------------------------------------------------------------- 00025 * Includes 00026 *----------------------------------------------------------------------------*/ 00027 00028 #include <string.h> 00029 #include <math.h> 00030 00031 #include <cpl.h> 00032 00033 #include "kmo_utils.h" 00034 #include "kmo_functions.h" 00035 #include "kmo_priv_reconstruct.h" 00036 #include "kmo_priv_functions.h" 00037 #include "kmo_priv_lcorr.h" 00038 #include "kmo_cpl_extensions.h" 00039 #include "kmo_dfs.h" 00040 #include "kmo_error.h" 00041 #include "kmo_utils.h" 00042 #include "kmo_constants.h" 00043 #include "kmo_debug.h" 00044 00045 /*----------------------------------------------------------------------------- 00046 * Functions prototypes 00047 *----------------------------------------------------------------------------*/ 00048 00049 static int kmos_reconstruct_check_inputs(cpl_frameset *, const char *); 00050 00051 static int kmos_reconstruct_create(cpl_plugin *); 00052 static int kmos_reconstruct_exec(cpl_plugin *); 00053 static int kmos_reconstruct_destroy(cpl_plugin *); 00054 static int kmos_reconstruct(cpl_parameterlist *, cpl_frameset *); 00055 00056 /*----------------------------------------------------------------------------- 00057 * Static variables 00058 *----------------------------------------------------------------------------*/ 00059 00060 static char kmos_reconstruct_description[] = 00061 "Data with or without noise is reconstructed into a cube using X/Y/LCAL, YCAL\n" 00062 "The input data can contain noise extensions and will be reconstructed into\n" 00063 "additional extensions.\n" 00064 "If an OH spectrum is given in the SOF file the lambda axis will be corrected\n" 00065 "using the OH lines as reference.\n" 00066 "\n" 00067 "---------------------------------------------------------------------------\n" 00068 " Input files:\n" 00069 "\n" 00070 " DO KMOS \n" 00071 " category Type Explanation Required #Frames\n" 00072 " -------- ----- ----------- -------- -------\n" 00073 " DARK or RAW/F2D data with Y 1 \n" 00074 " FLAT_ON or RAW/F2D or without noise \n" 00075 " ARC_ON or RAW/F2D \n" 00076 " OBJECT or RAW \n" 00077 " STD or RAW \n" 00078 " SCIENCE RAW \n" 00079 " XCAL F2D x-direction calib. frame Y 1 \n" 00080 " YCAL F2D y-direction calib. frame Y 1 \n" 00081 " LCAL F2D Wavelength calib. frame Y 1 \n" 00082 " WAVE_BAND F2L Table with start-/end-wavelengths Y 1 \n" 00083 " OH_SPEC F1S Vector holding OH lines N 1 \n" 00084 "\n" 00085 " Output files:\n" 00086 "\n" 00087 " DO KMOS\n" 00088 " category Type Explanation\n" 00089 " -------- ----- -----------\n" 00090 " CUBE_DARK or F3I Reconstructed cube \n" 00091 " CUBE_FLAT or RAW/F2D with or without noise\n" 00092 " CUBE_ARC or \n" 00093 " CUBE_OBJECT or \n" 00094 " CUBE_STD or \n" 00095 " CUBE_SCIENCE \n" 00096 "---------------------------------------------------------------------------\n" 00097 "\n"; 00098 00099 /*----------------------------------------------------------------------------- 00100 * Functions code 00101 *----------------------------------------------------------------------------*/ 00102 00103 /*----------------------------------------------------------------------------*/ 00107 /*----------------------------------------------------------------------------*/ 00108 00111 /*----------------------------------------------------------------------------*/ 00120 /*----------------------------------------------------------------------------*/ 00121 int cpl_plugin_get_info(cpl_pluginlist *list) 00122 { 00123 cpl_recipe *recipe = cpl_calloc(1, sizeof *recipe); 00124 cpl_plugin *plugin = &recipe->interface; 00125 00126 cpl_plugin_init(plugin, 00127 CPL_PLUGIN_API, 00128 KMOS_BINARY_VERSION, 00129 CPL_PLUGIN_TYPE_RECIPE, 00130 "kmos_reconstruct", 00131 "Performs the cube reconstruction", 00132 kmos_reconstruct_description, 00133 "Alex Agudo Berbel, Y. Jung", 00134 "usd-help@eso.org", 00135 kmos_get_license(), 00136 kmos_reconstruct_create, 00137 kmos_reconstruct_exec, 00138 kmos_reconstruct_destroy); 00139 00140 cpl_pluginlist_append(list, plugin); 00141 00142 return 0; 00143 } 00144 00145 /*----------------------------------------------------------------------------*/ 00153 /*----------------------------------------------------------------------------*/ 00154 static int kmos_reconstruct_create(cpl_plugin *plugin) 00155 { 00156 cpl_recipe *recipe; 00157 cpl_parameter *p; 00158 00159 /* Check that the plugin is part of a valid recipe */ 00160 if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE) 00161 recipe = (cpl_recipe *)plugin; 00162 else 00163 return -1; 00164 00165 /* Create the parameters list in the cpl_recipe object */ 00166 recipe->parameters = cpl_parameterlist_new(); 00167 00168 /* Fill the parameters list */ 00169 /* --imethod */ 00170 p = cpl_parameter_new_value("kmos.kmos_reconstruct.imethod", 00171 CPL_TYPE_STRING, 00172 "Method to use for interpolation. [\"NN\" (nearest neighbour), " 00173 "\"lwNN\" (linear weighted nearest neighbor), " 00174 "\"swNN\" (square weighted nearest neighbor), " 00175 "\"MS\" (Modified Shepard's method)" 00176 "\"CS\" (Cubic spline)]", 00177 "kmos.kmos_reconstruct", "CS"); 00178 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "imethod"); 00179 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV); 00180 cpl_parameterlist_append(recipe->parameters, p); 00181 00182 /* --neighborhoodRange */ 00183 p = cpl_parameter_new_value("kmos.kmos_reconstruct.neighborhoodRange", 00184 CPL_TYPE_DOUBLE, 00185 "Defines the range to search for neighbors. in pixels", 00186 "kmos.kmos_reconstruct", 1.001); 00187 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "neighborhoodRange"); 00188 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV); 00189 cpl_parameterlist_append(recipe->parameters, p); 00190 00191 /* --flux */ 00192 p = cpl_parameter_new_value("kmos.kmos_reconstruct.flux", CPL_TYPE_BOOL, 00193 "TRUE: Apply flux conservation. FALSE: otherwise", 00194 "kmos.kmos_reconstruct", FALSE); 00195 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "flux"); 00196 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV); 00197 cpl_parameterlist_append(recipe->parameters, p); 00198 00199 /* --detectorimage */ 00200 p = cpl_parameter_new_value("kmos.kmos_reconstruct.detectorimage", 00201 CPL_TYPE_BOOL, 00202 "TRUE: if resampled detector frame should be " 00203 "created, FALSE: otherwise", 00204 "kmos.kmos_reconstruct", FALSE); 00205 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "detimg"); 00206 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV); 00207 cpl_parameterlist_append(recipe->parameters, p); 00208 00209 /* --file_extension */ 00210 p = cpl_parameter_new_value("kmos.kmos_reconstruct.file_extension", 00211 CPL_TYPE_BOOL, 00212 "TRUE: if OBS_ID keyword should be appended to " 00213 "output frames, FALSE: otherwise", 00214 "kmos.kmos_reconstruct", FALSE); 00215 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "file_extension"); 00216 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV); 00217 cpl_parameterlist_append(recipe->parameters, p); 00218 00219 /* --pix_scale */ 00220 p = cpl_parameter_new_value("kmos.kmos_reconstruct.pix_scale", 00221 CPL_TYPE_DOUBLE, 00222 "Change the pixel scale [arcsec]. " 00223 "Default of 0.2\" results into cubes of 14x14pix, " 00224 "a scale of 0.1\" results into cubes of 28x28pix, etc.", 00225 "kmos.kmos_reconstruct", KMOS_PIX_RESOLUTION); 00226 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "pix_scale"); 00227 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV); 00228 cpl_parameterlist_append(recipe->parameters, p); 00229 00230 /* --xcal_interpolation */ 00231 p = cpl_parameter_new_value("kmos.kmos_reconstruct.xcal_interpolation", 00232 CPL_TYPE_BOOL, 00233 "TRUE: Interpolate xcal between rotator angles. FALSE: otherwise", 00234 "kmos.kmos_reconstruct", TRUE); 00235 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "xcal_interpolation"); 00236 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV); 00237 cpl_parameterlist_append(recipe->parameters, p); 00238 00239 /* Add parameters for band-definition */ 00240 kmos_band_pars_create(recipe->parameters, "kmos.kmos_reconstruct"); 00241 return 0; 00242 } 00243 00244 /*----------------------------------------------------------------------------*/ 00250 /*----------------------------------------------------------------------------*/ 00251 static int kmos_reconstruct_exec(cpl_plugin *plugin) 00252 { 00253 cpl_recipe *recipe; 00254 00255 /* Get the recipe out of the plugin */ 00256 if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE) 00257 recipe = (cpl_recipe *)plugin; 00258 else return -1; 00259 00260 return kmos_reconstruct(recipe->parameters, recipe->frames); 00261 } 00262 00263 /*----------------------------------------------------------------------------*/ 00269 /*----------------------------------------------------------------------------*/ 00270 static int kmos_reconstruct_destroy(cpl_plugin *plugin) 00271 { 00272 cpl_recipe *recipe; 00273 00274 /* Get the recipe out of the plugin */ 00275 if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE) 00276 recipe = (cpl_recipe *)plugin; 00277 else return -1 ; 00278 00279 cpl_parameterlist_delete(recipe->parameters); 00280 return 0 ; 00281 } 00282 00283 /*----------------------------------------------------------------------------*/ 00296 /*----------------------------------------------------------------------------*/ 00297 static int kmos_reconstruct(cpl_parameterlist *parlist, cpl_frameset *frameset) 00298 { 00299 const cpl_parameter * par ; 00300 const char * imethod ; 00301 int flux, detectorimage, xcal_interpolation, 00302 file_extension ; 00303 double neighborhoodRange, pix_scale, scaling ; 00304 cpl_frame * input_frame ; 00305 cpl_frame * xcal_frame ; 00306 cpl_frame * ycal_frame ; 00307 cpl_frame * lcal_frame ; 00308 cpl_frame * ref_spectrum_frame ; 00309 float * pdet_img_data ; 00310 float * pdet_img_noise ; 00311 float * slice ; 00312 const char * input_frame_name ; 00313 const char * output_frame_name ; 00314 const char * filter_id ; 00315 char * keyword ; 00316 char * suffix ; 00317 char * obs_suffix ; 00318 char * extname ; 00319 int * bounds ; 00320 cpl_image * det_img_data ; 00321 cpl_image * det_img_noise; 00322 cpl_image * tmp_img ; 00323 cpl_imagelist * cube_data ; 00324 cpl_imagelist * cube_noise ; 00325 cpl_propertylist * main_header ; 00326 cpl_propertylist * sub_header ; 00327 cpl_propertylist * tmp_header ; 00328 cpl_table * band_table ; 00329 gridDefinition gd ; 00330 cpl_polynomial * lcorr_coeffs ; 00331 main_fits_desc desc1 ; 00332 int i, j, index, ifu_nr, detImgCube, l, x, y ; 00333 00334 /* Initialise */ 00335 detImgCube = FALSE ; 00336 00337 /* Check Entries */ 00338 if (parlist == NULL || frameset == NULL) { 00339 cpl_msg_error(__func__, "Null Inputs") ; 00340 cpl_error_set(__func__, CPL_ERROR_NULL_INPUT) ; 00341 return -1 ; 00342 } 00343 00344 /* Get parameters */ 00345 par = cpl_parameterlist_find_const(parlist,"kmos.kmos_reconstruct.imethod"); 00346 imethod = cpl_parameter_get_string(par) ; 00347 par = cpl_parameterlist_find_const(parlist, "kmos.kmos_reconstruct.flux"); 00348 flux = cpl_parameter_get_bool(par); 00349 par = cpl_parameterlist_find_const(parlist, 00350 "kmos.kmos_reconstruct.detectorimage"); 00351 detectorimage = cpl_parameter_get_bool(par); 00352 par = cpl_parameterlist_find_const(parlist, 00353 "kmos.kmos_reconstruct.neighborhoodRange"); 00354 neighborhoodRange = cpl_parameter_get_double(par) ; 00355 kmos_band_pars_load(parlist, "kmos.kmos_reconstruct"); 00356 par = cpl_parameterlist_find_const(parlist, 00357 "kmos.kmos_reconstruct.file_extension"); 00358 file_extension = cpl_parameter_get_bool(par); 00359 par = cpl_parameterlist_find_const(parlist, 00360 "kmos.kmos_reconstruct.pix_scale"); 00361 pix_scale = cpl_parameter_get_double(par) ; 00362 par = cpl_parameterlist_find_const(parlist, 00363 "kmos.kmos_reconstruct.xcal_interpolation"); 00364 xcal_interpolation = cpl_parameter_get_bool(par); 00365 00366 /* Check Parameters */ 00367 if (strcmp(imethod, "NN") && strcmp(imethod, "lwNN") && 00368 strcmp(imethod, "swNN") && strcmp(imethod, "MS") && 00369 strcmp(imethod, "CS")) { 00370 cpl_msg_error(__func__, 00371 "imethod must be 'NN','lwNN','swNN','MS' or 'CS'") ; 00372 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00373 return -1 ; 00374 } 00375 if (neighborhoodRange <= 0.0) { 00376 cpl_msg_error(__func__, 00377 "neighborhoodRange must be greater than 0.0") ; 00378 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00379 return -1 ; 00380 } 00381 if (pix_scale < 0.01 || pix_scale > 0.4) { 00382 cpl_msg_error(__func__, "pix_scale must be between 0.01 and 0.4"); 00383 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00384 return -1 ; 00385 } 00386 00387 /* Identify the RAW and CALIB frames in the input frameset */ 00388 if (kmo_dfs_set_groups(frameset, "kmos_reconstruct") != 1) { 00389 cpl_msg_error(__func__, "Cannot identify RAW and CALIB frames") ; 00390 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00391 return -1 ; 00392 } 00393 00394 /* IO File names */ 00395 if (cpl_frameset_count_tags(frameset, DARK) == 1) { 00396 input_frame_name = DARK; 00397 output_frame_name = CUBE_DARK; 00398 } else if (cpl_frameset_count_tags(frameset, FLAT_ON) == 1) { 00399 input_frame_name = FLAT_ON; 00400 output_frame_name = CUBE_FLAT; 00401 } else if (cpl_frameset_count_tags(frameset, ARC_ON) == 1) { 00402 input_frame_name = ARC_ON; 00403 output_frame_name = CUBE_ARC; 00404 } else if (cpl_frameset_count_tags(frameset, OBJECT) == 1) { 00405 input_frame_name = OBJECT; 00406 output_frame_name = CUBE_OBJECT; 00407 } else if (cpl_frameset_count_tags(frameset, STD) == 1) { 00408 input_frame_name = STD; 00409 output_frame_name = CUBE_STD; 00410 } else if (cpl_frameset_count_tags(frameset, SCIENCE) == 1) { 00411 input_frame_name = SCIENCE; 00412 output_frame_name = CUBE_SCIENCE; 00413 } else { 00414 cpl_msg_error(__func__, "Missing Inputs") ; 00415 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00416 return -1 ; 00417 } 00418 00419 /* Check the inputs consistency */ 00420 if (kmos_reconstruct_check_inputs(frameset, input_frame_name) != 1) { 00421 cpl_msg_error(__func__, "Input frameset is not consistent") ; 00422 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00423 return -1 ; 00424 } 00425 00426 /* Instrument setup */ 00427 input_frame = kmo_dfs_get_frame(frameset, input_frame_name); 00428 suffix = kmo_dfs_get_suffix(input_frame, TRUE, TRUE); 00429 cpl_msg_info(__func__, "Detected instrument setup: %s", suffix+1); 00430 cpl_free(suffix); 00431 00432 /* Get frames */ 00433 xcal_frame = kmo_dfs_get_frame(frameset, XCAL); 00434 ycal_frame = kmo_dfs_get_frame(frameset, YCAL) ; 00435 lcal_frame = kmo_dfs_get_frame(frameset, LCAL) ; 00436 ref_spectrum_frame = kmo_dfs_get_frame(frameset, OH_SPEC) ; 00437 00438 /* Get input frame primary header */ 00439 main_header = kmo_dfs_load_primary_header(frameset, input_frame_name); 00440 00441 /* Use the OBS ID in the file name */ 00442 if (file_extension) { 00443 obs_suffix = cpl_sprintf("_%d", 00444 cpl_propertylist_get_int(main_header, OBS_ID)); 00445 } else { 00446 obs_suffix = cpl_sprintf("%s", ""); 00447 } 00448 00449 /* Save Product (reconstructed) Main Header */ 00450 kmo_dfs_save_main_header(frameset, output_frame_name, obs_suffix, 00451 input_frame, NULL, parlist, cpl_func); 00452 00453 /* Save Product (detector image) Main Header */ 00454 if (detectorimage == TRUE) 00455 kmo_dfs_save_main_header(frameset, DET_IMG_REC, obs_suffix, 00456 input_frame, NULL, parlist, cpl_func); 00457 00458 /* Grid definition, wavelength start and end are set later */ 00459 kmclipm_setup_grid(&gd, imethod, neighborhoodRange, pix_scale, 0.); 00460 00461 /* READ Wavelength bounds */ 00462 tmp_header = kmo_dfs_load_primary_header(frameset, XCAL); 00463 bounds = kmclipm_extract_bounds(tmp_header); 00464 cpl_propertylist_delete(tmp_header); 00465 00466 kmo_init_fits_desc(&desc1); 00467 desc1 = kmo_identify_fits_header(cpl_frame_get_filename(input_frame)); 00468 00469 /* Loop through detectors */ 00470 for (i = 1; i <= KMOS_NR_DETECTORS; i++) { 00471 cpl_msg_info("","Processing detector No. %d", i); 00472 00473 /* Load LCAL image close to ROTANGLE 0. assuming that the */ 00474 /* wavelength range doesn't differ with different ROTANGLEs. */ 00475 print_cal_angle_msg_once = FALSE; 00476 print_xcal_angle_msg_once = FALSE; 00477 if (i==1) { 00478 print_cal_angle_msg_once = TRUE; 00479 print_xcal_angle_msg_once = TRUE; 00480 } 00481 00482 /* Read Filter ID */ 00483 keyword = cpl_sprintf("%s%d%s", IFU_FILTID_PREFIX,i,IFU_FILTID_POSTFIX); 00484 filter_id = cpl_propertylist_get_string(main_header, keyword); 00485 cpl_free(keyword); 00486 00487 /* Load Band Info */ 00488 band_table = kmo_dfs_load_table(frameset, WAVE_BAND, 1, FALSE); 00489 kmclipm_setup_grid_band_lcal(&gd, filter_id, band_table); 00490 cpl_table_delete(band_table); 00491 00492 /* Create empty detector images */ 00493 if (detectorimage == TRUE) { 00494 det_img_data = cpl_image_new( 00495 gd.x.dim*gd.y.dim*KMOS_IFUS_PER_DETECTOR, gd.l.dim, 00496 CPL_TYPE_FLOAT); 00497 pdet_img_data = cpl_image_get_data_float(det_img_data); 00498 00499 det_img_noise = cpl_image_new( 00500 gd.x.dim*gd.y.dim*KMOS_IFUS_PER_DETECTOR, gd.l.dim, 00501 CPL_TYPE_FLOAT); 00502 pdet_img_noise = cpl_image_get_data_float(det_img_noise); 00503 } 00504 00505 /* Loop on IFUs */ 00506 for (j = 0; j < KMOS_IFUS_PER_DETECTOR; j++) { 00507 ifu_nr = (i-1)*KMOS_IFUS_PER_DETECTOR + j + 1; 00508 00509 /* Load sub-header*/ 00510 sub_header = kmo_dfs_load_sub_header(frameset, input_frame_name, 00511 i, FALSE); 00512 00513 /* Check if IFU is valid according to main header keywords */ 00514 if (bounds[2*(ifu_nr-1)] != -1 && bounds[2*(ifu_nr-1)+1] != -1) { 00515 /* IFU valid */ 00516 00517 /* calc WCS & update subheader */ 00518 kmo_calc_wcs_gd(main_header, sub_header, ifu_nr, gd); 00519 kmclipm_update_property_int(sub_header, NAXIS, 3, 00520 "number of data axes"); 00521 kmclipm_update_property_int(sub_header, NAXIS1, 00522 gd.x.dim, "length of data axis 1"); 00523 kmclipm_update_property_int(sub_header, NAXIS2, 00524 gd.y.dim, "length of data axis 2"); 00525 kmclipm_update_property_int(sub_header, NAXIS3, 00526 gd.l.dim, "length of data axis 3"); 00527 00528 /* Reconstruct data and noise (if available) */ 00529 kmo_reconstruct_sci(ifu_nr, bounds[2*(ifu_nr-1)], 00530 bounds[2*(ifu_nr-1)+1], input_frame, input_frame_name, 00531 NULL, NULL, NULL, xcal_frame, ycal_frame, lcal_frame, 00532 NULL, NULL, &gd, &cube_data, &cube_noise, flux, 00533 0, xcal_interpolation); 00534 00535 /* Reconstruct again using OH_SPEC */ 00536 if (ref_spectrum_frame != NULL && cube_data != NULL) { 00537 lcorr_coeffs = kmo_lcorr_get(cube_data, sub_header, 00538 ref_spectrum_frame, gd, filter_id, ifu_nr); 00539 cpl_imagelist_delete(cube_data); 00540 if (cube_noise != NULL) cpl_imagelist_delete(cube_noise); 00541 kmo_reconstruct_sci(ifu_nr, bounds[2*(ifu_nr-1)], 00542 bounds[2*(ifu_nr-1)+1], input_frame, 00543 input_frame_name, NULL, NULL, NULL, xcal_frame, 00544 ycal_frame, lcal_frame, lcorr_coeffs, NULL, &gd, 00545 &cube_data, &cube_noise, flux, 0, 00546 xcal_interpolation); 00547 cpl_polynomial_delete(lcorr_coeffs); 00548 } 00549 00550 /* Scale flux according to pixel_scale */ 00551 tmp_img = cpl_imagelist_get(cube_data, 0); 00552 scaling = (cpl_image_get_size_x(tmp_img)* 00553 cpl_image_get_size_y(tmp_img)) / 00554 (KMOS_SLITLET_X*KMOS_SLITLET_Y); 00555 cpl_imagelist_divide_scalar(cube_data, scaling); 00556 if (cube_noise != NULL) 00557 cpl_imagelist_divide_scalar(cube_noise, scaling); 00558 } else { 00559 /* IFU invalid */ 00560 cube_data = cube_noise = NULL ; 00561 } 00562 00563 /* Fill detector images */ 00564 if (detectorimage) { 00565 if (cube_data != NULL) { 00566 for (l = 0; l < gd.l.dim; l++) { 00567 slice = cpl_image_get_data_float( 00568 cpl_imagelist_get(cube_data, l)); 00569 for (y = 0; y < gd.y.dim; y++) { 00570 for (x = 0; x < gd.x.dim; x++) { 00571 int ix = x + 00572 y* gd.x.dim + 00573 j* gd.x.dim*gd.y.dim + //IFU offset 00574 l* gd.x.dim*gd.y.dim*KMOS_IFUS_PER_DETECTOR; 00575 pdet_img_data[ix] = slice[x + y*gd.x.dim]; 00576 } 00577 } 00578 } 00579 } 00580 if (cube_noise != NULL) { 00581 detImgCube = TRUE; 00582 for (l = 0; l < gd.l.dim; l++) { 00583 slice = cpl_image_get_data_float( 00584 cpl_imagelist_get(cube_noise, l)); 00585 for (y = 0; y < gd.y.dim; y++) { 00586 for (x = 0; x < gd.x.dim; x++) { 00587 int ix = x + 00588 y* gd.x.dim + 00589 j* gd.x.dim*gd.y.dim + //IFU offset 00590 l* gd.x.dim*gd.y.dim*KMOS_IFUS_PER_DETECTOR; 00591 pdet_img_noise[ix] = slice[x + y*gd.x.dim]; 00592 } 00593 } 00594 } 00595 } 00596 } 00597 00598 /* Save Product (reconstructed) Data Cube */ 00599 extname = kmo_extname_creator(ifu_frame, ifu_nr, EXT_DATA); 00600 kmclipm_update_property_string(sub_header, EXTNAME, extname, 00601 "FITS extension name"); 00602 cpl_free(extname); 00603 kmo_dfs_save_cube(cube_data, output_frame_name, obs_suffix, 00604 sub_header, 0./0.); 00605 00606 /* Save Product (reconstructed) Noise Cube */ 00607 if (cube_noise != NULL) { 00608 extname = kmo_extname_creator(ifu_frame, ifu_nr, EXT_NOISE); 00609 kmclipm_update_property_string(sub_header, EXTNAME, 00610 extname, "FITS extension name"); 00611 cpl_free(extname); 00612 kmo_dfs_save_cube(cube_noise, output_frame_name, obs_suffix, 00613 sub_header, 0./0.); 00614 } 00615 00616 cpl_imagelist_delete(cube_data); 00617 if (cube_noise != NULL) cpl_imagelist_delete(cube_noise); 00618 cpl_propertylist_delete(sub_header); 00619 } // IFUs loop 00620 00621 if (detectorimage) { 00622 /* Save Product (detector image) Data Cube */ 00623 index = kmo_identify_index(cpl_frame_get_filename(input_frame), i, 00624 FALSE); 00625 tmp_header = kmclipm_propertylist_load( 00626 cpl_frame_get_filename(input_frame), index); 00627 kmo_save_det_img_ext(det_img_data, gd, i, DET_IMG_REC, 00628 obs_suffix, tmp_header, FALSE, FALSE); 00629 cpl_propertylist_delete(tmp_header); 00630 cpl_image_delete(det_img_data); 00631 00632 if (detImgCube) { 00633 /* Save Product (detector image) Noise Cube */ 00634 /* Index changes (*2) if Noise extensions are there */ 00635 if (desc1.ex_noise) { 00636 index = kmo_identify_index( 00637 cpl_frame_get_filename(input_frame), i, TRUE); 00638 } 00639 tmp_header = kmclipm_propertylist_load( 00640 cpl_frame_get_filename(input_frame), index); 00641 kmo_save_det_img_ext(det_img_noise, gd, i, DET_IMG_REC, 00642 obs_suffix, tmp_header, FALSE, TRUE); 00643 cpl_propertylist_delete(tmp_header); 00644 } 00645 cpl_image_delete(det_img_noise); 00646 } 00647 } // Detectors loop 00648 cpl_propertylist_delete(main_header); 00649 cpl_free(obs_suffix); 00650 kmo_free_fits_desc(&desc1); 00651 cpl_free(bounds); 00652 00653 return 0; 00654 } 00655 00658 /*----------------------------------------------------------------------------*/ 00665 /*----------------------------------------------------------------------------*/ 00666 static int kmos_reconstruct_check_inputs( 00667 cpl_frameset * frameset, 00668 const char * input_frame_name) 00669 { 00670 cpl_propertylist * lcal_header ; 00671 cpl_propertylist * input_header ; 00672 int nb_xcal, nb_ycal, nb_lcal, nb_oh_spec ; 00673 char * keyword ; 00674 const char * filter_id ; 00675 const char * filter_id_tmp ; 00676 int i ; 00677 00678 /* Check entries */ 00679 if (frameset == NULL || input_frame_name == NULL) return -1; 00680 00681 /* Count frames */ 00682 if (cpl_frameset_count_tags(frameset, DARK) != 1 && 00683 cpl_frameset_count_tags(frameset, FLAT_ON) != 1 && 00684 cpl_frameset_count_tags(frameset, ARC_ON) != 1 && 00685 cpl_frameset_count_tags(frameset, OBJECT) != 1 && 00686 cpl_frameset_count_tags(frameset, STD) != 1 && 00687 cpl_frameset_count_tags(frameset, SCIENCE) != 1) { 00688 cpl_msg_error(__func__, "1 data frame must be provided") ; 00689 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00690 return 0 ; 00691 } 00692 nb_xcal = cpl_frameset_count_tags(frameset, XCAL) ; 00693 nb_ycal = cpl_frameset_count_tags(frameset, YCAL) ; 00694 nb_lcal = cpl_frameset_count_tags(frameset, LCAL) ; 00695 if (nb_xcal != 1 || nb_ycal != 1 || nb_lcal != 1) { 00696 cpl_msg_error(__func__, "Exactly 1 XCAL/YCAL/LCAL expected") ; 00697 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00698 return 0 ; 00699 } 00700 if (cpl_frameset_count_tags(frameset, WAVE_BAND) != 1) { 00701 cpl_msg_error(__func__, "Missing WAVE_BAND") ; 00702 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00703 return 0 ; 00704 } 00705 nb_oh_spec = cpl_frameset_count_tags(frameset, OH_SPEC) ; 00706 if (nb_oh_spec != 0 && nb_oh_spec != 1) { 00707 cpl_msg_error(__func__, "Only 1 reference spectrum can be provided") ; 00708 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00709 return 0 ; 00710 } 00711 00712 /* Check filters, rotation consistencies */ 00713 kmo_check_frame_setup(frameset, XCAL, YCAL, TRUE, FALSE, TRUE); 00714 kmo_check_frame_setup(frameset, XCAL, LCAL, TRUE, FALSE, TRUE); 00715 if (nb_oh_spec == 1) kmo_check_oh_spec_setup(frameset, XCAL); 00716 00717 if (cpl_frameset_count_tags(frameset, DARK) != 1) { 00718 kmo_check_frame_setup(frameset, XCAL, input_frame_name, TRUE, FALSE, 00719 FALSE); 00720 } 00721 kmo_check_frame_setup_md5_xycal(frameset); 00722 kmo_check_frame_setup_md5(frameset); 00723 00724 /* Check filters */ 00725 input_header = kmo_dfs_load_primary_header(frameset, input_frame_name); 00726 lcal_header = kmo_dfs_load_primary_header(frameset, LCAL); 00727 for (i = 1; i <= KMOS_NR_DETECTORS; i++) { 00728 // ESO INS FILTi ID 00729 keyword = cpl_sprintf("%s%d%s", IFU_FILTID_PREFIX,i,IFU_FILTID_POSTFIX); 00730 filter_id = cpl_propertylist_get_string(lcal_header, keyword); 00731 if (strcmp(filter_id, "IZ") && strcmp(filter_id, "YJ") && 00732 strcmp(filter_id, "H") && strcmp(filter_id, "K") && 00733 strcmp(filter_id, "HK")) { 00734 cpl_free(keyword) ; 00735 cpl_propertylist_delete(input_header) ; 00736 cpl_propertylist_delete(lcal_header) ; 00737 cpl_msg_error(__func__,"LCAL Filter ID must be IZ, YJ, H, HK or K"); 00738 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00739 return 0 ; 00740 } 00741 00742 /* Dark not taken with filter */ 00743 if (!strcmp(input_frame_name, DARK)) { 00744 filter_id_tmp = cpl_propertylist_get_string(input_header, keyword); 00745 if (strcmp(filter_id, filter_id_tmp)) { 00746 cpl_free(keyword) ; 00747 cpl_propertylist_delete(input_header) ; 00748 cpl_propertylist_delete(lcal_header) ; 00749 cpl_msg_error(__func__,"Filter mismatch"); 00750 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ; 00751 return 0 ; 00752 } 00753 } 00754 cpl_free(keyword); 00755 } 00756 cpl_propertylist_delete(input_header) ; 00757 cpl_propertylist_delete(lcal_header) ; 00758 00759 return 1 ; 00760 } 00761