GRAVI Pipeline Reference Manual  0.7.11
gravi_dual.c
1 /* $Id: gravi_dual.c,v 1.29 2011/12/3 09:16:12 nazouaoui Exp $
2  *
3  * This file is part of the GRAVI Pipeline
4  * Copyright (C) 2002,2003 European Southern Observatory
5  *
6  * This program is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 2 of the License, or
9  * (at your option) any later version.
10  *
11  * This program is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with this program; if not, write to the Free Software
18  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19  */
20 
21 /*
22  * $Author: nazouaoui $
23  * $Date: 2011/12/3 09:16:12 $
24  * $Revision: 1.29 $
25  * $Name: $
26  */
27 
28 #ifdef HAVE_CONFIG_H
29 #include <config.h>
30 #endif
31 
32 /*-----------------------------------------------------------------------------
33  Includes
34  -----------------------------------------------------------------------------*/
35 
36 #include <cpl.h>
37 #include <stdio.h>
38 #include <string.h>
39 
40 #include "gravi_utils.h"
41 #include "gravi_pfits.h"
42 #include "gravi_dfs.h"
43 #include "gravi_calib.h"
44 #include "gravi_vis.h"
45 #include "gravi_data.h"
46 
47 
48 
49 /*-----------------------------------------------------------------------------
50  Private function prototypes
51  -----------------------------------------------------------------------------*/
52 
53 static int gravi_dual_create(cpl_plugin *);
54 static int gravi_dual_exec(cpl_plugin *);
55 static int gravi_dual_destroy(cpl_plugin *);
56 static int gravi_dual(cpl_frameset *, const cpl_parameterlist *);
57 
58 /*-----------------------------------------------------------------------------
59  Static variables
60  -----------------------------------------------------------------------------*/
61 
62 static char gravi_dual_description[] =
63 "This recipe is associated to the observing template in dual field mode.\n"
64 "Its aim is to reduce the raw interferometric data acquired on calibrator or \n"
65 "science targets and to compute the intrumental visibility.\n"
66 "The output OIFITS file contains the values of the complex visibilities, \n"
67 "squared visibilities and the closure phases.\n"
68 "Description DO category\n"
69 "Required input :\n"
70 " Master dark " DARK "\n"
71 " Master flat " FLAT "\n"
72 " Bad pixel map " BAD "\n"
73 " Master wave " WAVE "\n"
74 " Master p2vm " P2VM "\n"
75 " Raw observation file " GRAVI_DUAL_CALIB " or "
76 GRAVI_DUAL_SCIENCE "\n"
77 "Optional input\n"
78 " Sky file (one or more) " DUAL_CAL_SKY_RAW " or "
79 DUAL_SCI_SKY_RAW "\n"
80 
81 "Output :\n"
82 " Master sky (one per input sky) " SKY "\n"
83 " Raw Visibilities " VIS_DUAL_CALIB " or " VIS_DUAL_SCIENCE "\n"
84 "Optional output\n"
85 " Preprocessed files (--preproc-file) " PREPROC "\n"
86 " Non averaged visibilities (--p2vmreduced-file) " P2VMREDUCED "\n"
87 "\n";
88 
89 /*-----------------------------------------------------------------------------
90  Function code
91  -----------------------------------------------------------------------------*/
92 
93 /*----------------------------------------------------------------------------*/
103 /*----------------------------------------------------------------------------*/
104 int cpl_plugin_get_info(cpl_pluginlist * list)
105 {
106  cpl_recipe * recipe = cpl_calloc(1, sizeof *recipe );
107  cpl_plugin * plugin = &recipe->interface;
108 
109  if (cpl_plugin_init(plugin,
110  CPL_PLUGIN_API,
111  GRAVI_BINARY_VERSION,
112  CPL_PLUGIN_TYPE_RECIPE,
113  "gravi_dual",
114  "This recipe is used to compute the squared visibility",
115  gravi_dual_description,
116  "Firstname Lastname",
117  PACKAGE_BUGREPORT,
118  gravi_get_license(),
119  gravi_dual_create,
120  gravi_dual_exec,
121  gravi_dual_destroy)) {
122  cpl_msg_error(cpl_func, "Plugin initialization failed");
123  (void)cpl_error_set_where(cpl_func);
124  return 1;
125  }
126 
127  if (cpl_pluginlist_append(list, plugin)) {
128  cpl_msg_error(cpl_func, "Error adding plugin to list");
129  (void)cpl_error_set_where(cpl_func);
130  return 1;
131  }
132 
133  return 0;
134 }
135 
136 /*----------------------------------------------------------------------------*/
144 /*----------------------------------------------------------------------------*/
145 static int gravi_dual_create(cpl_plugin * plugin)
146 {
147  cpl_recipe * recipe;
148  cpl_parameter * p;
149 
150  /* Do not create the recipe if an error code is already set */
151  if (cpl_error_get_code() != CPL_ERROR_NONE) {
152  cpl_msg_error(cpl_func, "%s():%d: An error is already set: %s",
153  cpl_func, __LINE__, cpl_error_get_where());
154  return (int)cpl_error_get_code();
155  }
156 
157  if (plugin == NULL) {
158  cpl_msg_error(cpl_func, "Null plugin");
159  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
160  }
161 
162  /* Verify plugin type */
163  if (cpl_plugin_get_type(plugin) != CPL_PLUGIN_TYPE_RECIPE) {
164  cpl_msg_error(cpl_func, "Plugin is not a recipe");
165  cpl_ensure_code(0, (int)CPL_ERROR_TYPE_MISMATCH);
166  }
167 
168  /* Get the recipe */
169  recipe = (cpl_recipe *)plugin;
170 
171  /* Create the parameters list in the cpl_recipe object */
172  recipe->parameters = cpl_parameterlist_new();
173  if (recipe->parameters == NULL) {
174  cpl_msg_error(cpl_func, "Parameter list allocation failed");
175  cpl_ensure_code(0, (int)CPL_ERROR_ILLEGAL_OUTPUT);
176  }
177 
178  /* Fill the parameters list */
179 
180  /* --Save the preproc files */
181  p = cpl_parameter_new_value("gravi.preproc_param.preproc_file",
182  CPL_TYPE_BOOL, "Save the preprocessed file", "gravi.gravi_all_p2vm", FALSE);
183  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "preproc-file");
184  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
185  cpl_parameterlist_append(recipe->parameters, p);
186 
187  /* --Save the p2vm reduced files */
188  p = cpl_parameter_new_value("gravi.preproc_param.p2vmreduced_file",
189  CPL_TYPE_BOOL, "Save the file with non averaged visibilities", "gravi.gravi_dual", FALSE);
190  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "p2vmreduced-file");
191  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
192  cpl_parameterlist_append(recipe->parameters, p);
193 
194  /* --Save the preproc files */
195  p = cpl_parameter_new_value("gravi.vis_reduce.vis_correction",
196  CPL_TYPE_BOOL, "Correction of visibilities looses", "gravi.gravi_dual", FALSE);
197  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "vis-correction");
198  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
199  cpl_parameterlist_append(recipe->parameters, p);
200 
201  /* --SNR threshold for fringe DET in FT */
202  p = cpl_parameter_new_value("gravi.snr_threshold_ft",
203  CPL_TYPE_DOUBLE, "SNR threshold for fringe detection in FT (>0)", "gravi", 3.0);
204  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "snr-threshold-ft");
205  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
206  cpl_parameterlist_append(recipe->parameters, p);
207 
208  /* --Minimum locking ratio to accept SC frame */
209  p = cpl_parameter_new_value("gravi.fringedet_threshold_sc",
210  CPL_TYPE_DOUBLE, "Fringe-detection ratio threshold to accept SC frame (0..1)", "gravi", 0.8);
211  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "fringedet-threshold-sc");
212  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
213  cpl_parameterlist_append(recipe->parameters, p);
214 
215  /* --Debias the V2 of SC */
216  p = cpl_parameter_new_value("gravi.debias_sc",
217  CPL_TYPE_BOOL, "Subtrast the V2 bias from SC", "gravi", TRUE);
218  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "debias-sc");
219  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
220  cpl_parameterlist_append(recipe->parameters, p);
221 
222  /* --Debias the V2 of FT */
223  p = cpl_parameter_new_value("gravi.debias_ft",
224  CPL_TYPE_BOOL, "Subtrast the V2 bias from FT", "gravi", TRUE);
225  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "debias-ft");
226  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
227  cpl_parameterlist_append(recipe->parameters, p);
228 
229  /* --Use the metrology */
230  p = cpl_parameter_new_value("gravi.use_met",
231  CPL_TYPE_BOOL, "Use the metrology to rephase the SC data in real-time", "gravi", FALSE);
232  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "use-met");
233  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
234  cpl_parameterlist_append(recipe->parameters, p);
235 
236  /* --Use the vis_flat */
237  p = cpl_parameter_new_value("gravi.vis_flat",
238  CPL_TYPE_BOOL, "Flat the OIVIS by the VIS_FLAT file", "gravi", TRUE);
239  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "vis-flat");
240  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
241  cpl_parameterlist_append(recipe->parameters, p);
242 
243  return 0;
244 }
245 
246 /*----------------------------------------------------------------------------*/
252 /*----------------------------------------------------------------------------*/
253 static int gravi_dual_exec(cpl_plugin * plugin)
254 {
255 
256  cpl_recipe * recipe;
257  int recipe_status;
258  cpl_errorstate initial_errorstate = cpl_errorstate_get();
259 
260 
261  /* Return immediately if an error code is already set */
262  if (cpl_error_get_code() != CPL_ERROR_NONE) {
263  cpl_msg_error(cpl_func, "%s():%d: An error is already set: %s",
264  cpl_func, __LINE__, cpl_error_get_where());
265  return (int)cpl_error_get_code();
266  }
267 
268  if (plugin == NULL) {
269  cpl_msg_error(cpl_func, "Null plugin");
270  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
271  }
272 
273  /* Verify plugin type */
274  if (cpl_plugin_get_type(plugin) != CPL_PLUGIN_TYPE_RECIPE) {
275  cpl_msg_error(cpl_func, "Plugin is not a recipe");
276  cpl_ensure_code(0, (int)CPL_ERROR_TYPE_MISMATCH);
277  }
278 
279  /* Get the recipe */
280  recipe = (cpl_recipe *)plugin;
281 
282  /* Verify parameter and frame lists */
283  if (recipe->parameters == NULL) {
284  cpl_msg_error(cpl_func, "Recipe invoked with NULL parameter list");
285  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
286  }
287  if (recipe->frames == NULL) {
288  cpl_msg_error(cpl_func, "Recipe invoked with NULL frame set");
289  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
290  }
291 
292  /* Invoke the recipe */
293  recipe_status = gravi_dual(recipe->frames, recipe->parameters);
294 
295  /* Ensure DFS-compliance of the products */
296 
297  if (cpl_dfs_update_product_header(recipe->frames)) {
298  if (!recipe_status){
299  recipe_status = (int)cpl_error_get_code();
300  }
301  }
302 
303  if (!cpl_errorstate_is_equal(initial_errorstate)) {
304  /* Dump the error history since recipe execution start.
305  At this point the recipe cannot recover from the error */
306  cpl_errorstate_dump(initial_errorstate, CPL_FALSE, NULL);
307  }
308 
309  return recipe_status;
310 }
311 
312 /*----------------------------------------------------------------------------*/
318 /*----------------------------------------------------------------------------*/
319 static int gravi_dual_destroy(cpl_plugin * plugin)
320 {
321  cpl_recipe * recipe;
322 
323  if (plugin == NULL) {
324  cpl_msg_error(cpl_func, "Null plugin");
325  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
326  }
327 
328  /* Verify plugin type */
329  if (cpl_plugin_get_type(plugin) != CPL_PLUGIN_TYPE_RECIPE) {
330  cpl_msg_error(cpl_func, "Plugin is not a recipe");
331  cpl_ensure_code(0, (int)CPL_ERROR_TYPE_MISMATCH);
332  }
333 
334  /* Get the recipe */
335  recipe = (cpl_recipe *)plugin;
336 
337  cpl_parameterlist_delete(recipe->parameters);
338 
339  return 0;
340 }
341 
342 
343 /*----------------------------------------------------------------------------*/
351 /*----------------------------------------------------------------------------*/
352 static int gravi_dual(cpl_frameset * frameset,
353  const cpl_parameterlist * parlist)
354 {
355  cpl_frameset * recipe_frameset=NULL, * wave_frameset=NULL, * dark_frameset=NULL,
356  * darkcalib_frameset=NULL, * skyframeset=NULL,
357  * gain_frameset=NULL, * p2vm_frameset=NULL, * badpix_frameset=NULL, *used_frameset=NULL,
358  * visflat_frameset=NULL;
359  cpl_propertylist * applist=NULL, * primary_hdr=NULL, * primary_hdr_dark=NULL, * wave_plist=NULL
360  , *bad_primary_hdr=NULL;
361  cpl_table * detector_table=NULL;
362  cpl_frame * frame=NULL;
363  const cpl_parameter * p=NULL, * p_=NULL;
364  cpl_parameter * par=NULL;
365  char * out_dark=NULL, * preproc_name=NULL;
366  const char * frame_tag=NULL, * filename=NULL, * temp=NULL;
367  gravi_data * p2vm_map=NULL, * data=NULL, * darkOrSky_map=NULL, * wave_map=NULL, * dark_map=NULL,
368  * profile_map=NULL, * badpix_map=NULL, * vis_flat=NULL;
369  gravi_data ** calib_datas=NULL;
370  gravi_data * preproc_data=NULL;
371  gravi_data ** raw_data=NULL, * p2vm_reduce=NULL;
372  gravi_data * oi_vis=NULL;
373  darkOrSky_map = NULL, wave_map = NULL;
374  int nb_frame, i, j, nwave_ft, nwave_sc, n_region;
375  int * shutter=NULL;
376  int nb_vis;
377  int testDarkSky = 0;
378  char * vis_name=NULL;
379 
380  /* Message */
381  cpl_msg_set_time_on();
382  cpl_msg_set_component_on();
383  cpl_msg_info(cpl_func,"Start function -- cleanup");
384 
385 
386  /*
387  * Identify the RAW and CALIB frames in the input frameset
388  */
389 
390  cpl_ensure_code(gravi_dfs_set_groups(frameset) == CPL_ERROR_NONE,
391  cpl_error_get_code()) ;
392 
393  p2vm_frameset = gravi_frameset_extract_p2vm_file(frameset);
394 
395  /* - Extract a set of frame p2vm and set of frame dark wave */
396  recipe_frameset = gravi_frameset_extract_dual_data(frameset);
397  darkcalib_frameset = gravi_frameset_extract_dark_file(frameset);
398  wave_frameset = gravi_frameset_extract_wave_file(frameset);
399  dark_frameset = gravi_frameset_extract_dark(frameset);
400  gain_frameset = gravi_frameset_extract_flat_file(frameset);
401  badpix_frameset = gravi_frameset_extract_badpix(frameset);
402  skyframeset = gravi_frameset_extract_dual_sky_data (frameset);
403  visflat_frameset = gravi_frameset_extract_visflat(frameset);
404 
405  /* To use this recipe the frameset must contain at least
406  * one p2vm frame. */
407  if (cpl_frameset_is_empty(recipe_frameset)) {
408  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE,
409  "No DUAL frame on the frameset");
410  goto cleanup;
411  }
412 
413  /* To use this recipe the frameset must contain the p2vm, wave and
414  * gain calibration file. */
415  if ((cpl_frameset_is_empty(p2vm_frameset)) ||
416  (cpl_frameset_is_empty(wave_frameset)) ||
417  (cpl_frameset_is_empty(gain_frameset))) {
418  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE,
419  "No P2VM or FLAT or WAVE file on the frameset") ;
420  goto cleanup;
421  }
422 
423  /* Get the number of the p2vm frame contained in the frameset */
424  nb_frame = cpl_frameset_get_size(recipe_frameset);
425 
426  /*
427  * Identify the DARK in the input frameset
428  */
429 
430  if (!cpl_frameset_is_empty(dark_frameset)) {
431 
432  frame = cpl_frameset_get_position(dark_frameset, 0);
433  filename = cpl_frame_get_filename(frame);
434 
435  cpl_msg_info (cpl_func, "File %s is a raw DARK file", FILESHORT(filename));
436  data = gravi_data_load(filename);
437 
438  /* The dark file must contains all the shutter close */
439  primary_hdr = gravi_data_get_propertylist(data, GRAVI_PRIMARY_HDR_NAME_EXT);
440  shutter = gravi_shutters_check(primary_hdr);
441  if (!((shutter[0] == 0) && (shutter[1] == 0) &&
442  (shutter[2] == 0) && (shutter[3] == 0))){
443  cpl_error_set_message(cpl_func, CPL_ERROR_ILLEGAL_OUTPUT, "Shutter problem in the DARK");
444  goto cleanup;
445  }
446 
447  /* Compute dark */
448  dark_map = gravi_compute_dark(data);
449 
450  CPLCHECK_CLEAN("Could not compute the DARK map");
451 
452  /* Save the dark map */
453  primary_hdr_dark = gravi_data_get_propertylist(dark_map,
454  GRAVI_PRIMARY_HDR_NAME_EXT);
455  applist = gravi_propertylist_get_qc (primary_hdr_dark);
456  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, DARK);
457 
458  // out_dark = cpl_sprintf("master_dark.%s", FILESHORT(filename));
459  out_dark = gravi_data_product_name (filename, "dark");
460 
461  gravi_data_save(dark_map, frameset, out_dark, parlist,
462  dark_frameset, frame, "gravi_dual", applist);
463 
464  CPLCHECK_CLEAN("Could not save the DARK map");
465 
466  FREE (cpl_propertylist_delete,applist);
467  FREE (gravi_data_delete,data);
468  FREE (cpl_free,shutter);
469 
470  testDarkSky = 1;
471  }
472  else if (!cpl_frameset_is_empty(darkcalib_frameset)) {
473 
474  frame = cpl_frameset_get_position(darkcalib_frameset, 0);
475  filename = cpl_frame_get_filename(frame);
476 
477  cpl_msg_info (cpl_func, "File %s is a master DARK file already computed", FILESHORT(filename));
478  dark_map = gravi_data_load(filename);
479  cpl_frameset_insert (dark_frameset, cpl_frame_duplicate (frame));
480 
481  testDarkSky = 1;
482  }
483  else
484  cpl_msg_info (cpl_func, "There is no DARK in the frame set");
485 
486  /*
487  * Identify the SKY in the input frameset
488  */
489 
490  if (!cpl_frameset_is_empty(skyframeset)) {
491 
492  /* Load the raw SKY */
493  frame = cpl_frameset_get_position(skyframeset, 0);
494  filename = cpl_frame_get_filename(frame);
495 
496  cpl_msg_info (cpl_func, "File %s is a raw SKY file", FILESHORT(filename));
497 
498  data = gravi_data_load(filename);
499  primary_hdr = gravi_data_get_propertylist(data, GRAVI_PRIMARY_HDR_NAME_EXT);
500 
501  /* Compute the SKY map */
502  darkOrSky_map = gravi_compute_dark(data);
503  FREE (gravi_data_delete,data);
504 
505  CPLCHECK_CLEAN("Error while computing the darkOrSky_map");
506 
507  /* Save the SKY map */
508  cpl_frameset_insert (dark_frameset, cpl_frame_duplicate (frame));
509 
510  primary_hdr_dark = gravi_data_get_propertylist(darkOrSky_map,
511  GRAVI_PRIMARY_HDR_NAME_EXT);
512 
513  applist = gravi_propertylist_get_qc (primary_hdr_dark);
514 
515  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, SKY);
516 
517  out_dark = gravi_data_product_name (filename, "sky");
518 
519  gravi_data_save(darkOrSky_map, frameset, out_dark, parlist,
520  skyframeset, frame, "gravi_dual", applist);
521 
522  CPLCHECK_CLEAN("Could not save the darkOrSky_map");
523 
524  FREE (cpl_propertylist_delete,applist);
525  FREE (cpl_free,out_dark);
526  }
527  else if (testDarkSky) {
528  darkOrSky_map = dark_map;
529  dark_map = NULL;
530  testDarkSky = 0;
531  cpl_msg_info (cpl_func, "There is no sky in the frame set");
532  }
533  else {
534  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE, "No dark or sky frames in the frameset");
535  goto cleanup;
536  }
537 
538  FREE (cpl_frameset_delete,skyframeset);
539 
540  /*
541  * Identify the BAD in the input frameset
542  */
543 
544  if (!cpl_frameset_is_empty(badpix_frameset)){
545 
546  /* Load the BAD pixel map */
547  frame = cpl_frameset_get_position(badpix_frameset, 0);
548  filename = cpl_frame_get_filename(frame);
549 
550  cpl_msg_info (cpl_func, "File %s is a BAD pixel map", FILESHORT(filename));
551  badpix_map = gravi_data_load(filename);
552  }
553  else {
554  /* Compute the BAD pixel map from the DARK or SKY */
555  badpix_map = gravi_compute_badpix(darkOrSky_map, parlist);
556 
557  CPLCHECK_CLEAN("Could not compute the BAD pixel map");
558 
559  /* Save the BAD pixel map */
560  bad_primary_hdr = gravi_data_get_propertylist (badpix_map, GRAVI_PRIMARY_HDR_NAME_EXT);
561  applist = gravi_propertylist_get_qc (bad_primary_hdr);
562 
563  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, BAD);
564 
565  frame = cpl_frameset_get_position (dark_frameset, 0);
566 
567  gravi_data_save (badpix_map, frameset, "gravi_bad_map.fits", parlist,
568  dark_frameset, frame, "gravi_dual", applist);
569 
570  CPLCHECK_CLEAN("Could not save the BAD pixel map");
571 
572  FREE (cpl_propertylist_delete, applist);
573  }
574 
575  /*
576  * Identify the FLAT in the input frameset
577  */
578 
579  frame = cpl_frameset_get_position(gain_frameset, 0);
580  filename = cpl_frame_get_filename(frame);
581 
582  cpl_msg_info (cpl_func, "File %s is a master FLAT already computed", FILESHORT(filename));
583  profile_map = gravi_data_load(filename);
584 
585  /*
586  * Identify the WAVE in the input frameset
587  */
588 
589  frame = cpl_frameset_get_position(wave_frameset, 0);
590  filename = cpl_frame_get_filename(frame);
591 
592  cpl_msg_info (cpl_func, "File %s is a master WAVE map already computed", FILESHORT(filename));
593  wave_map = gravi_data_load(filename);
594 
595  /* Check if exist in the input frameset */
596  if ((wave_map == NULL) || (profile_map == NULL)) {
597  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE, "No WAVE or FLAT frames in the frameset");
598  goto cleanup;
599  }
600 
601  /*
602  * Identify the P2VM in the input frameset
603  */
604 
605  frame = cpl_frameset_get_position(p2vm_frameset, 0);
606  filename = cpl_frame_get_filename(frame);
607 
608  cpl_msg_info (cpl_func, "File %s is a P2VM file already computed", FILESHORT(filename));
609  p2vm_map = gravi_data_load (filename);
610 
611  CPLCHECK_CLEAN("Error while loading the P2VM map");
612 
613  /*
614  * Identify the VIS_FLAT in the input frameset
615  */
616  if (!cpl_frameset_is_empty(visflat_frameset))
617  {
618  frame = cpl_frameset_get_position(visflat_frameset, 0);
619  filename = cpl_frame_get_filename(frame);
620  cpl_msg_info (cpl_func, "File %s is a VIS_FLAT already computed", FILESHORT(filename));
621  vis_flat = gravi_data_load(filename);
622  }
623 
624  /*
625  * Prepare all the input and calibration for the reduction
626  */
627 
628  /* Check if the data with the p2vm have the same spectral n */
629  detector_table = gravi_data_get_table (p2vm_map, GRAVI_IMAGING_DETECTOR_FT_EXT);
630  n_region = cpl_table_get_nrow(detector_table);
631  if (n_region > 24) {
632  wave_plist = gravi_data_get_oi_propertylist
633  (p2vm_map, GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_SC),
634  GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_SC,0,2));
635  nwave_sc = gravi_pfits_get_nwave (wave_plist);
636  wave_plist = gravi_data_get_oi_propertylist (p2vm_map,
637  GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_FT),
638  GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_FT,0,2));
639  nwave_ft = gravi_pfits_get_nwave (wave_plist);
640  }
641  else {
642  wave_plist = gravi_data_get_oi_propertylist
643  (p2vm_map, GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_SC), GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_SC,0,1));
644  nwave_sc = gravi_pfits_get_nwave (wave_plist);
645  wave_plist = gravi_data_get_oi_propertylist (p2vm_map,
646  GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_FT), GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_FT,0,1));
647  nwave_ft = gravi_pfits_get_nwave (wave_plist);
648  }
649 
650  primary_hdr = gravi_data_get_propertylist (wave_map, GRAVI_PRIMARY_HDR_NAME_EXT);
651 
652  CPLCHECK_CLEAN("Missed one the preproc parameters");
653 
654  /* Construction of the calibration data */
655  int nb_calfile;
656  if (testDarkSky)
657  nb_calfile = 5;
658  else
659  nb_calfile = 4;
660 
661  calib_datas = cpl_malloc(nb_calfile * sizeof(gravi_data *));
662 
663  calib_datas[0] = darkOrSky_map;
664  calib_datas[1] = wave_map;
665  calib_datas[2] = profile_map;
666  calib_datas[3] = badpix_map;
667  if (testDarkSky)
668  calib_datas[4] = dark_map;
669 
670  /* insert calibration frame into the used frameset */
671  used_frameset=cpl_frameset_new();
672  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(dark_frameset));
673  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(badpix_frameset));
674  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(wave_frameset));
675  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(gain_frameset));
676  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(p2vm_frameset));
677 
678  p = cpl_parameterlist_find_const(parlist, "gravi.preproc_param.preproc_file");
679 
680  p_ = cpl_parameterlist_find_const(parlist, "gravi.preproc_param.p2vmreduced_file");
681 
682  /* Construction of the preproc data and extract the spectra using
683  * the results of the spectral calibration */
684 
685  nb_vis = 0;
686  oi_vis = gravi_data_new(0);
687  primary_hdr = gravi_data_get_propertylist (p2vm_map,
688  GRAVI_PRIMARY_HDR_NAME_EXT);
689 
690  const char * resol_p2vm = gravi_pfits_get_resolution (primary_hdr);
691  const char * resol_file;
692 
693 
694  /*
695  * Loop on input RAW frames to be reduced
696  */
697 
698  int thread;
699  for (thread = 0; thread < nb_frame; thread++){
700 
701  frame = cpl_frameset_get_position (recipe_frameset, thread);
702  frame_tag = cpl_frame_get_tag (frame);
703  filename = cpl_frame_get_filename (frame);
704 
705  data = gravi_data_load(filename);
706  primary_hdr = gravi_data_get_propertylist (data,
707  GRAVI_PRIMARY_HDR_NAME_EXT);
708 
709  /* Check the shutter for extracting the files can be used to compute
710  * the visibilities */
711  shutter = gravi_shutters_check(primary_hdr);
712  resol_file = gravi_pfits_get_resolution (primary_hdr);
713 
714  /* Check shutters */
715  if ((shutter [0] != 1) || (shutter [1] != 1) ||
716  (shutter [2] != 1) || (shutter [3] != 1)) {
717  cpl_msg_warning( cpl_func, "RAW data shutters: %d-%d-%d-%d",
718  shutter[0], shutter[1], shutter[2], shutter[3] );
719  }
720 
721  cpl_msg_info (cpl_func, "Preprocessing of the dual"
722  " field file %s", FILESHORT(filename));
723 
724  preproc_data = gravi_preproc(data, calib_datas, nb_calfile, p2vm_map, parlist);
725 
726  FREE (gravi_data_delete,data);
727  FREE (cpl_free,shutter);
728 
729  CPLCHECK_CLEAN("Cannot preproc the data");
730 
731  /* Option save the proproc file */
732  if (cpl_parameter_get_bool(p)){
733  // preproc_name = cpl_sprintf("gravi_preproc.%s", FILESHORT(filename));
734  preproc_name = gravi_data_product_name (filename, "preproc");
735  cpl_frameset * preproc_frame = cpl_frameset_new ();
736  cpl_frameset_insert (preproc_frame, cpl_frame_duplicate (frame));
737  cpl_frameset_join(preproc_frame, cpl_frameset_duplicate(used_frameset)); // add the calibration frameset
738 
739  applist = gravi_propertylist_get_qc (gravi_data_get_propertylist
740  (preproc_data, GRAVI_PRIMARY_HDR_NAME_EXT));
741 
742  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, PREPROC);
743 
744  gravi_data_save(preproc_data, frameset, preproc_name, parlist,
745  preproc_frame, frame, "gravi_dual", applist);
746 
747  CPLCHECK_CLEAN("Cannot save the PREPROC product");
748 
749  FREE (cpl_free,preproc_name);
750  FREE (cpl_propertylist_delete,applist);
751  FREE (cpl_frameset_delete,preproc_frame);
752  }
753 
754  /* Compute the flux and visibilities for each telescope and
755  * per acquisition with the P2VM applied to preproc_data */
756  cpl_msg_info (cpl_func, "Application of the p2vm to file %s", FILESHORT(filename));
757  p2vm_reduce = gravi_p2vm_reduce (preproc_data, p2vm_map, "gravi_dual", parlist);
758  FREE (gravi_data_delete,preproc_data);
759 
760  CPLCHECK_CLEAN("Cannot apply p2vm to the preproc data");
761 
762  /*
763  * Compute the QC0 about tau0 from piezo signals
764  */
765  gravi_compute_tau0 (p2vm_reduce);
766 
767  /* Visibility and flux are averaged and the followings
768  * are saved in Visibility data in tables VIS, VIS2 and T3 */
769  cpl_msg_info (cpl_func, "Computing the average visibilities of file %s", FILESHORT(filename));
770  gravi_vis_reduce (oi_vis, p2vm_reduce, nb_vis, parlist, "gravi_dual");
771 
772  CPLCHECK_CLEAN("Cannot average the visibilities");
773 
774  /* Compute the QC parameters of the TF
775  * TODO: to be done only if a calibration star */
776  cpl_msg_info( cpl_func, "Compute the QC TF parameters" );
777  cpl_errorstate prestate = cpl_errorstate_get();
778 
779  gravi_data * oi_tf = gravi_compute_tf( oi_vis );
780 
781  cpl_propertylist_copy_property_regexp( gravi_data_get_propertylist (oi_vis, GRAVI_PRIMARY_HDR_NAME_EXT),
782  gravi_data_get_propertylist (oi_tf, GRAVI_PRIMARY_HDR_NAME_EXT),
783  ".*QC TF.*", 0);
784  FREE (gravi_data_delete,oi_tf);
785 
786  /* If an error is catch when computing the QC parameters, dump this error but continue */
787  if( cpl_error_get_code() ) {
788  cpl_msg_warning(cpl_func, "Cannot compute the QC TF parameters for this observation");
789  cpl_errorstate_dump (prestate, 0, NULL);
790  cpl_error_reset();
791  }
792 
793  /* Save the p2vmreduced file */
794  if (cpl_parameter_get_bool(p_)){
795  cpl_msg_info( cpl_func, "Save the p2vmreduced file" );
796 
797  // char * p2vm_name = cpl_sprintf("dual_reduced.%s", FILESHORT(filename));
798  char * p2vm_name = gravi_data_product_name (filename, "p2vmreduced");
799 
800  cpl_frameset * p2vmreduced_frame = cpl_frameset_new ();
801  cpl_frameset_insert (p2vmreduced_frame, cpl_frame_duplicate (frame));
802  cpl_frameset_join(p2vmreduced_frame, cpl_frameset_duplicate(used_frameset));
803  applist = gravi_propertylist_get_qc (gravi_data_get_propertylist
804  (p2vm_reduce, GRAVI_PRIMARY_HDR_NAME_EXT));
805 
806  if ((strcmp(frame_tag, DUAL_CAL_SKY_RAW) == 0) ||
807  (strcmp(frame_tag, DUAL_SCI_SKY_RAW) == 0))
808  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, SKYREDUCED);
809  else
810  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, P2VMREDUCED);
811 
812  gravi_data_save(p2vm_reduce, frameset, p2vm_name, parlist,
813  p2vmreduced_frame, frame, "gravi_dual", applist);
814 
815  CPLCHECK_CLEAN("Cannot save the P2VMREDUCED product");
816 
817  FREE (cpl_free,p2vm_name);
818  FREE (cpl_propertylist_delete,applist);
819  FREE (cpl_frameset_delete,p2vmreduced_frame);
820  }
821 
822  nb_vis++;
823  cpl_msg_info(cpl_func,"Free the p2vmreduced");
824  FREE (gravi_data_delete,p2vm_reduce);
825  }
826  /* End loop on the input files to reduce */
827 
828  /* Eventually performn the flat with VIS_FLAT if provided */
829  if (gravi_param_get_bool(parlist, "gravi.vis_flat")) {
830  if (vis_flat)
831  {
832  cpl_msg_info (cpl_func, "Perform the visibility flat with the provided VIS_FLAT file");
833  gravi_flat_vis (oi_vis, vis_flat);
834  CPLCHECK_CLEAN("Cannot flat the OI_VIS");
835  }
836  else
837  {
838  cpl_msg_warning (cpl_func,"vis-flat option requested but no VIS_FLAT in the frameset");
839  }
840  }
841 
842  /* Save the output data file */
843  frame = cpl_frameset_get_position (recipe_frameset, 0);
844  frame_tag = cpl_frame_get_tag (frame);
845  primary_hdr = gravi_data_get_propertylist(oi_vis,
846  GRAVI_PRIMARY_HDR_NAME_EXT);
847  applist = gravi_propertylist_get_qc (primary_hdr);
848 
849  cpl_propertylist_append_int(applist, GRAVI_NIGHT_OBS,
850  cpl_propertylist_get_int (primary_hdr, GRAVI_NIGHT_OBS));
851  cpl_propertylist_append_string(applist, "DATE-OBS",
852  cpl_propertylist_get_string (primary_hdr, "DATE-OBS"));
853 
854 
855  if ((strcmp(frame_tag, GRAVI_DUAL_CALIB) == 0)) {
856  vis_name = gravi_data_product_name (filename, "dualcal");
857  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_DUAL_CALIB);
858  }
859  else if ((strcmp(frame_tag, GRAVI_DUAL_SCIENCE) == 0)) {
860  vis_name = gravi_data_product_name (filename, "dualsci");
861  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_DUAL_SCIENCE);
862  }
863  else {
864  vis_name = gravi_data_product_name (filename, "dualsci");
865  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_DUAL_SCIENCE);
866  }
867 
868  cpl_frameset_join (recipe_frameset, cpl_frameset_duplicate (used_frameset));
869 
870  gravi_data_save(oi_vis, frameset, vis_name, parlist,
871  recipe_frameset, frame, "gravi_dual", applist);
872 
873  CPLCHECK_CLEAN("Cannot save the VIS product");
874 
875 
876  /* Terminate the function */
877  goto cleanup;
878 
879 cleanup:
880  /* Deallocation of all variables */
881  cpl_msg_info(cpl_func,"Memory cleanup");
882 
883  FREE (gravi_data_delete,dark_map);
884  FREE (gravi_data_delete,data);
885  FREE (cpl_free,shutter);
886  FREE (cpl_free,preproc_name);
887  FREE (gravi_data_delete,preproc_data);
888  FREE (gravi_data_delete,profile_map);
889  FREE (gravi_data_delete,darkOrSky_map);
890  FREE (gravi_data_delete,wave_map);
891  FREE (gravi_data_delete,badpix_map);
892  FREE (gravi_data_delete,p2vm_map);
893  FREE (gravi_data_delete,vis_flat);
894  FREE (cpl_frameset_delete,darkcalib_frameset);
895  FREE (cpl_frameset_delete,wave_frameset);
896  FREE (cpl_frameset_delete,dark_frameset);
897  FREE (cpl_frameset_delete,gain_frameset);
898  FREE (cpl_frameset_delete,badpix_frameset);
899  FREE (gravi_data_delete,p2vm_reduce);
900  FREE (cpl_free,calib_datas);
901  FREE (cpl_free,vis_name);
902  FREE (cpl_propertylist_delete,applist);
903  FREE (gravi_data_delete,oi_vis);
904  FREE (cpl_frameset_delete,p2vm_frameset);
905  FREE (cpl_frameset_delete,visflat_frameset);
906  FREE (cpl_frameset_delete,recipe_frameset);
907  FREE (cpl_frameset_delete,used_frameset);
908 
909  cpl_msg_info(cpl_func,"Exit function");
910  return (int)cpl_error_get_code();
911 }
912 
913