GRAVI Pipeline Reference Manual  0.7.12
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  /* --Correction of visibility losses */
195  p = cpl_parameter_new_value("gravi.vis_reduce.vis_correction",
196  CPL_TYPE_BOOL, "Correction of visibility losses", "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  /* --STATE threshold for fringe DET in FT */
209  p = cpl_parameter_new_value("gravi.state_threshold_ft",
210  CPL_TYPE_DOUBLE, "STATE threshold for fringe detection in FT (>=0)", "gravi.vis_reduce", 1.0);
211  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "state-threshold-ft");
212  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
213  cpl_parameterlist_append(recipe->parameters, p);
214 
215  /* --Minimum detection ratio to accept SC frame */
216  p = cpl_parameter_new_value("gravi.fringedet_threshold_sc",
217  CPL_TYPE_DOUBLE, "Fringe-detection ratio threshold to accept SC frame (0..1)", "gravi", 0.8);
218  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "fringedet-threshold-sc");
219  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
220  cpl_parameterlist_append(recipe->parameters, p);
221 
222  /* --vFactor threshold to accept SC frame */
223  p = cpl_parameter_new_value("gravi.vfactor_threshold_sc",
224  CPL_TYPE_DOUBLE, "vFactor threshold to accept SC frame (0..1)", "gravi.vis_reduce", 0.1);
225  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "vfactor-threshold-sc");
226  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
227  cpl_parameterlist_append(recipe->parameters, p);
228 
229  /* --Debias the V2 of SC */
230  p = cpl_parameter_new_value("gravi.debias_sc",
231  CPL_TYPE_BOOL, "Subtract the V2 bias from SC", "gravi", TRUE);
232  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "debias-sc");
233  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
234  cpl_parameterlist_append(recipe->parameters, p);
235 
236  /* --Debias the V2 of FT */
237  p = cpl_parameter_new_value("gravi.debias_ft",
238  CPL_TYPE_BOOL, "Subtract the V2 bias from FT", "gravi", TRUE);
239  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "debias-ft");
240  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
241  cpl_parameterlist_append(recipe->parameters, p);
242 
243  /* --Use the metrology */
244  p = cpl_parameter_new_value("gravi.use_met",
245  CPL_TYPE_BOOL, "Use the metrology to rephase the SC data in real-time", "gravi", FALSE);
246  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "use-met");
247  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
248  cpl_parameterlist_append(recipe->parameters, p);
249 
250  /* --Use the vis_flat */
251  p = cpl_parameter_new_value("gravi.vis_flat",
252  CPL_TYPE_BOOL, "Flat the OIVIS by the VIS_FLAT file", "gravi", TRUE);
253  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "vis-flat");
254  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
255  cpl_parameterlist_append(recipe->parameters, p);
256 
257  return 0;
258 }
259 
260 /*----------------------------------------------------------------------------*/
266 /*----------------------------------------------------------------------------*/
267 static int gravi_dual_exec(cpl_plugin * plugin)
268 {
269 
270  cpl_recipe * recipe;
271  int recipe_status;
272  cpl_errorstate initial_errorstate = cpl_errorstate_get();
273 
274 
275  /* Return immediately if an error code is already set */
276  if (cpl_error_get_code() != CPL_ERROR_NONE) {
277  cpl_msg_error(cpl_func, "%s():%d: An error is already set: %s",
278  cpl_func, __LINE__, cpl_error_get_where());
279  return (int)cpl_error_get_code();
280  }
281 
282  if (plugin == NULL) {
283  cpl_msg_error(cpl_func, "Null plugin");
284  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
285  }
286 
287  /* Verify plugin type */
288  if (cpl_plugin_get_type(plugin) != CPL_PLUGIN_TYPE_RECIPE) {
289  cpl_msg_error(cpl_func, "Plugin is not a recipe");
290  cpl_ensure_code(0, (int)CPL_ERROR_TYPE_MISMATCH);
291  }
292 
293  /* Get the recipe */
294  recipe = (cpl_recipe *)plugin;
295 
296  /* Verify parameter and frame lists */
297  if (recipe->parameters == NULL) {
298  cpl_msg_error(cpl_func, "Recipe invoked with NULL parameter list");
299  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
300  }
301  if (recipe->frames == NULL) {
302  cpl_msg_error(cpl_func, "Recipe invoked with NULL frame set");
303  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
304  }
305 
306  /* Invoke the recipe */
307  recipe_status = gravi_dual(recipe->frames, recipe->parameters);
308 
309  /* Ensure DFS-compliance of the products */
310 
311  if (cpl_dfs_update_product_header(recipe->frames)) {
312  if (!recipe_status){
313  recipe_status = (int)cpl_error_get_code();
314  }
315  }
316 
317  if (!cpl_errorstate_is_equal(initial_errorstate)) {
318  /* Dump the error history since recipe execution start.
319  At this point the recipe cannot recover from the error */
320  cpl_errorstate_dump(initial_errorstate, CPL_FALSE, NULL);
321  }
322 
323  return recipe_status;
324 }
325 
326 /*----------------------------------------------------------------------------*/
332 /*----------------------------------------------------------------------------*/
333 static int gravi_dual_destroy(cpl_plugin * plugin)
334 {
335  cpl_recipe * recipe;
336 
337  if (plugin == NULL) {
338  cpl_msg_error(cpl_func, "Null plugin");
339  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
340  }
341 
342  /* Verify plugin type */
343  if (cpl_plugin_get_type(plugin) != CPL_PLUGIN_TYPE_RECIPE) {
344  cpl_msg_error(cpl_func, "Plugin is not a recipe");
345  cpl_ensure_code(0, (int)CPL_ERROR_TYPE_MISMATCH);
346  }
347 
348  /* Get the recipe */
349  recipe = (cpl_recipe *)plugin;
350 
351  cpl_parameterlist_delete(recipe->parameters);
352 
353  return 0;
354 }
355 
356 
357 /*----------------------------------------------------------------------------*/
365 /*----------------------------------------------------------------------------*/
366 static int gravi_dual(cpl_frameset * frameset,
367  const cpl_parameterlist * parlist)
368 {
369  cpl_frameset * recipe_frameset=NULL, * wave_frameset=NULL, * dark_frameset=NULL,
370  * darkcalib_frameset=NULL, * skyframeset=NULL,
371  * gain_frameset=NULL, * p2vm_frameset=NULL, * badpix_frameset=NULL, *used_frameset=NULL,
372  * visflat_frameset=NULL;
373  cpl_propertylist * applist=NULL, * primary_hdr=NULL, * primary_hdr_dark=NULL, * wave_plist=NULL
374  , *bad_primary_hdr=NULL;
375  cpl_table * detector_table=NULL;
376  cpl_frame * frame=NULL;
377  const cpl_parameter * p=NULL, * p_=NULL;
378  cpl_parameter * par=NULL;
379  char * out_dark=NULL, * preproc_name=NULL;
380  const char * frame_tag=NULL, * filename=NULL, * temp=NULL;
381  gravi_data * p2vm_map=NULL, * data=NULL, * darkOrSky_map=NULL, * wave_map=NULL, * dark_map=NULL,
382  * profile_map=NULL, * badpix_map=NULL, * vis_flat=NULL;
383  gravi_data ** calib_datas=NULL;
384  gravi_data * preproc_data=NULL;
385  gravi_data ** raw_data=NULL, * p2vm_reduce=NULL;
386  gravi_data * oi_vis=NULL;
387  darkOrSky_map = NULL, wave_map = NULL;
388  int nb_frame, i, j, nwave_ft, nwave_sc, n_region;
389  int * shutter=NULL;
390  int nb_vis;
391  int testDarkSky = 0;
392  char * vis_name=NULL;
393 
394  /* Message */
395  cpl_msg_set_time_on();
396  cpl_msg_set_component_on();
397  cpl_msg_info(cpl_func,"Start function -- cleanup");
398 
399 
400  /*
401  * Identify the RAW and CALIB frames in the input frameset
402  */
403 
404  cpl_ensure_code(gravi_dfs_set_groups(frameset) == CPL_ERROR_NONE,
405  cpl_error_get_code()) ;
406 
407  p2vm_frameset = gravi_frameset_extract_p2vm_file(frameset);
408 
409  /* - Extract a set of frame p2vm and set of frame dark wave */
410  recipe_frameset = gravi_frameset_extract_dual_data(frameset);
411  darkcalib_frameset = gravi_frameset_extract_dark_file(frameset);
412  wave_frameset = gravi_frameset_extract_wave_file(frameset);
413  dark_frameset = gravi_frameset_extract_dark(frameset);
414  gain_frameset = gravi_frameset_extract_flat_file(frameset);
415  badpix_frameset = gravi_frameset_extract_badpix(frameset);
416  skyframeset = gravi_frameset_extract_dual_sky_data (frameset);
417  visflat_frameset = gravi_frameset_extract_visflat(frameset);
418 
419  /* To use this recipe the frameset must contain at least
420  * one p2vm frame. */
421  if (cpl_frameset_is_empty(recipe_frameset)) {
422  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE,
423  "No DUAL frame on the frameset");
424  goto cleanup;
425  }
426 
427  /* To use this recipe the frameset must contain the p2vm, wave and
428  * gain calibration file. */
429  if ((cpl_frameset_is_empty(p2vm_frameset)) ||
430  (cpl_frameset_is_empty(wave_frameset)) ||
431  (cpl_frameset_is_empty(gain_frameset))) {
432  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE,
433  "No P2VM or FLAT or WAVE file on the frameset") ;
434  goto cleanup;
435  }
436 
437  /* Get the number of the p2vm frame contained in the frameset */
438  nb_frame = cpl_frameset_get_size(recipe_frameset);
439 
440  /*
441  * Identify the DARK in the input frameset
442  */
443 
444  if (!cpl_frameset_is_empty(dark_frameset)) {
445 
446  frame = cpl_frameset_get_position(dark_frameset, 0);
447  filename = cpl_frame_get_filename(frame);
448 
449  cpl_msg_info (cpl_func, "File %s is a raw DARK file", FILESHORT(filename));
450  data = gravi_data_load(filename);
451 
452  /* The dark file must contains all the shutter close */
453  primary_hdr = gravi_data_get_propertylist(data, GRAVI_PRIMARY_HDR_NAME_EXT);
454  shutter = gravi_shutters_check(primary_hdr);
455  if (!((shutter[0] == 0) && (shutter[1] == 0) &&
456  (shutter[2] == 0) && (shutter[3] == 0))){
457  cpl_error_set_message(cpl_func, CPL_ERROR_ILLEGAL_OUTPUT, "Shutter problem in the DARK");
458  goto cleanup;
459  }
460 
461  /* Compute dark */
462  dark_map = gravi_compute_dark(data);
463 
464  CPLCHECK_CLEAN("Could not compute the DARK map");
465 
466  /* Save the dark map */
467  primary_hdr_dark = gravi_data_get_propertylist(dark_map,
468  GRAVI_PRIMARY_HDR_NAME_EXT);
469  applist = gravi_propertylist_get_qc (primary_hdr_dark);
470  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, DARK);
471 
472  // out_dark = cpl_sprintf("master_dark.%s", FILESHORT(filename));
473  out_dark = gravi_data_product_name (filename, "dark");
474 
475  gravi_data_save(dark_map, frameset, out_dark, parlist,
476  dark_frameset, frame, "gravi_dual", applist);
477 
478  CPLCHECK_CLEAN("Could not save the DARK map");
479 
480  FREE (cpl_propertylist_delete,applist);
481  FREE (gravi_data_delete,data);
482  FREE (cpl_free,shutter);
483 
484  testDarkSky = 1;
485  }
486  else if (!cpl_frameset_is_empty(darkcalib_frameset)) {
487 
488  frame = cpl_frameset_get_position(darkcalib_frameset, 0);
489  filename = cpl_frame_get_filename(frame);
490 
491  cpl_msg_info (cpl_func, "File %s is a master DARK file already computed", FILESHORT(filename));
492  dark_map = gravi_data_load(filename);
493  cpl_frameset_insert (dark_frameset, cpl_frame_duplicate (frame));
494 
495  testDarkSky = 1;
496  }
497  else
498  cpl_msg_info (cpl_func, "There is no DARK in the frame set");
499 
500  /*
501  * Identify the SKY in the input frameset
502  */
503 
504  if (!cpl_frameset_is_empty(skyframeset)) {
505 
506  /* Load the raw SKY */
507  frame = cpl_frameset_get_position(skyframeset, 0);
508  filename = cpl_frame_get_filename(frame);
509 
510  cpl_msg_info (cpl_func, "File %s is a raw SKY file", FILESHORT(filename));
511 
512  data = gravi_data_load(filename);
513  primary_hdr = gravi_data_get_propertylist(data, GRAVI_PRIMARY_HDR_NAME_EXT);
514 
515  /* Compute the SKY map */
516  darkOrSky_map = gravi_compute_dark(data);
517  FREE (gravi_data_delete,data);
518 
519  CPLCHECK_CLEAN("Error while computing the darkOrSky_map");
520 
521  /* Save the SKY map */
522  cpl_frameset_insert (dark_frameset, cpl_frame_duplicate (frame));
523 
524  primary_hdr_dark = gravi_data_get_propertylist(darkOrSky_map,
525  GRAVI_PRIMARY_HDR_NAME_EXT);
526 
527  applist = gravi_propertylist_get_qc (primary_hdr_dark);
528 
529  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, SKY);
530 
531  out_dark = gravi_data_product_name (filename, "sky");
532 
533  gravi_data_save(darkOrSky_map, frameset, out_dark, parlist,
534  skyframeset, frame, "gravi_dual", applist);
535 
536  CPLCHECK_CLEAN("Could not save the darkOrSky_map");
537 
538  FREE (cpl_propertylist_delete,applist);
539  FREE (cpl_free,out_dark);
540  }
541  else if (testDarkSky) {
542  darkOrSky_map = dark_map;
543  dark_map = NULL;
544  testDarkSky = 0;
545  cpl_msg_info (cpl_func, "There is no sky in the frame set");
546  }
547  else {
548  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE, "No dark or sky frames in the frameset");
549  goto cleanup;
550  }
551 
552  FREE (cpl_frameset_delete,skyframeset);
553 
554  /*
555  * Identify the BAD in the input frameset
556  */
557 
558  if (!cpl_frameset_is_empty(badpix_frameset)){
559 
560  /* Load the BAD pixel map */
561  frame = cpl_frameset_get_position(badpix_frameset, 0);
562  filename = cpl_frame_get_filename(frame);
563 
564  cpl_msg_info (cpl_func, "File %s is a BAD pixel map", FILESHORT(filename));
565  badpix_map = gravi_data_load(filename);
566  }
567  else {
568  /* Compute the BAD pixel map from the DARK or SKY */
569  badpix_map = gravi_compute_badpix(darkOrSky_map, parlist);
570 
571  CPLCHECK_CLEAN("Could not compute the BAD pixel map");
572 
573  /* Save the BAD pixel map */
574  bad_primary_hdr = gravi_data_get_propertylist (badpix_map, GRAVI_PRIMARY_HDR_NAME_EXT);
575  applist = gravi_propertylist_get_qc (bad_primary_hdr);
576 
577  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, BAD);
578 
579  frame = cpl_frameset_get_position (dark_frameset, 0);
580 
581  gravi_data_save (badpix_map, frameset, "gravi_bad_map.fits", parlist,
582  dark_frameset, frame, "gravi_dual", applist);
583 
584  CPLCHECK_CLEAN("Could not save the BAD pixel map");
585 
586  FREE (cpl_propertylist_delete, applist);
587  }
588 
589  /*
590  * Identify the FLAT in the input frameset
591  */
592 
593  frame = cpl_frameset_get_position(gain_frameset, 0);
594  filename = cpl_frame_get_filename(frame);
595 
596  cpl_msg_info (cpl_func, "File %s is a master FLAT already computed", FILESHORT(filename));
597  profile_map = gravi_data_load(filename);
598 
599  /*
600  * Identify the WAVE in the input frameset
601  */
602 
603  frame = cpl_frameset_get_position(wave_frameset, 0);
604  filename = cpl_frame_get_filename(frame);
605 
606  cpl_msg_info (cpl_func, "File %s is a master WAVE map already computed", FILESHORT(filename));
607  wave_map = gravi_data_load(filename);
608 
609  /* Check if exist in the input frameset */
610  if ((wave_map == NULL) || (profile_map == NULL)) {
611  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE, "No WAVE or FLAT frames in the frameset");
612  goto cleanup;
613  }
614 
615  /*
616  * Identify the P2VM in the input frameset
617  */
618 
619  frame = cpl_frameset_get_position(p2vm_frameset, 0);
620  filename = cpl_frame_get_filename(frame);
621 
622  cpl_msg_info (cpl_func, "File %s is a P2VM file already computed", FILESHORT(filename));
623  p2vm_map = gravi_data_load (filename);
624 
625  CPLCHECK_CLEAN("Error while loading the P2VM map");
626 
627  /*
628  * Identify the VIS_FLAT in the input frameset
629  */
630  if (!cpl_frameset_is_empty(visflat_frameset))
631  {
632  frame = cpl_frameset_get_position(visflat_frameset, 0);
633  filename = cpl_frame_get_filename(frame);
634  cpl_msg_info (cpl_func, "File %s is a VIS_FLAT already computed", FILESHORT(filename));
635  vis_flat = gravi_data_load(filename);
636  }
637 
638  /*
639  * Prepare all the input and calibration for the reduction
640  */
641 
642  /* Check if the data with the p2vm have the same spectral n */
643  detector_table = gravi_data_get_table (p2vm_map, GRAVI_IMAGING_DETECTOR_FT_EXT);
644  n_region = cpl_table_get_nrow(detector_table);
645  if (n_region > 24) {
646  wave_plist = gravi_data_get_oi_propertylist
647  (p2vm_map, GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_SC),
648  GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_SC,0,2));
649  nwave_sc = gravi_pfits_get_nwave (wave_plist);
650  wave_plist = gravi_data_get_oi_propertylist (p2vm_map,
651  GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_FT),
652  GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_FT,0,2));
653  nwave_ft = gravi_pfits_get_nwave (wave_plist);
654  }
655  else {
656  wave_plist = gravi_data_get_oi_propertylist
657  (p2vm_map, GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_SC), GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_SC,0,1));
658  nwave_sc = gravi_pfits_get_nwave (wave_plist);
659  wave_plist = gravi_data_get_oi_propertylist (p2vm_map,
660  GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_FT), GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_FT,0,1));
661  nwave_ft = gravi_pfits_get_nwave (wave_plist);
662  }
663 
664  primary_hdr = gravi_data_get_propertylist (wave_map, GRAVI_PRIMARY_HDR_NAME_EXT);
665 
666  CPLCHECK_CLEAN("Missed one the preproc parameters");
667 
668  /* Construction of the calibration data */
669  int nb_calfile;
670  if (testDarkSky)
671  nb_calfile = 5;
672  else
673  nb_calfile = 4;
674 
675  calib_datas = cpl_malloc(nb_calfile * sizeof(gravi_data *));
676 
677  calib_datas[0] = darkOrSky_map;
678  calib_datas[1] = wave_map;
679  calib_datas[2] = profile_map;
680  calib_datas[3] = badpix_map;
681  if (testDarkSky)
682  calib_datas[4] = dark_map;
683 
684  /* insert calibration frame into the used frameset */
685  used_frameset=cpl_frameset_new();
686  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(dark_frameset));
687  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(badpix_frameset));
688  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(wave_frameset));
689  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(gain_frameset));
690  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(p2vm_frameset));
691 
692  p = cpl_parameterlist_find_const(parlist, "gravi.preproc_param.preproc_file");
693 
694  p_ = cpl_parameterlist_find_const(parlist, "gravi.preproc_param.p2vmreduced_file");
695 
696  /* Construction of the preproc data and extract the spectra using
697  * the results of the spectral calibration */
698 
699  nb_vis = 0;
700  oi_vis = gravi_data_new(0);
701  primary_hdr = gravi_data_get_propertylist (p2vm_map,
702  GRAVI_PRIMARY_HDR_NAME_EXT);
703 
704  const char * resol_p2vm = gravi_pfits_get_resolution (primary_hdr);
705  const char * resol_file;
706 
707 
708  /*
709  * Loop on input RAW frames to be reduced
710  */
711 
712  int thread;
713  for (thread = 0; thread < nb_frame; thread++){
714 
715  frame = cpl_frameset_get_position (recipe_frameset, thread);
716  frame_tag = cpl_frame_get_tag (frame);
717  filename = cpl_frame_get_filename (frame);
718 
719  data = gravi_data_load(filename);
720  primary_hdr = gravi_data_get_propertylist (data,
721  GRAVI_PRIMARY_HDR_NAME_EXT);
722 
723  /* Check the shutter for extracting the files can be used to compute
724  * the visibilities */
725  shutter = gravi_shutters_check(primary_hdr);
726  resol_file = gravi_pfits_get_resolution (primary_hdr);
727 
728  /* Check shutters */
729  if ((shutter [0] != 1) || (shutter [1] != 1) ||
730  (shutter [2] != 1) || (shutter [3] != 1)) {
731  cpl_msg_warning( cpl_func, "RAW data shutters: %d-%d-%d-%d",
732  shutter[0], shutter[1], shutter[2], shutter[3] );
733  }
734 
735  cpl_msg_info (cpl_func, "Preprocessing of the dual"
736  " field file %s", FILESHORT(filename));
737 
738  preproc_data = gravi_preproc(data, calib_datas, nb_calfile, p2vm_map, parlist);
739 
740  FREE (gravi_data_delete,data);
741  FREE (cpl_free,shutter);
742 
743  CPLCHECK_CLEAN("Cannot preproc the data");
744 
745  /* Option save the proproc file */
746  if (cpl_parameter_get_bool(p)){
747  // preproc_name = cpl_sprintf("gravi_preproc.%s", FILESHORT(filename));
748  preproc_name = gravi_data_product_name (filename, "preproc");
749  cpl_frameset * preproc_frame = cpl_frameset_new ();
750  cpl_frameset_insert (preproc_frame, cpl_frame_duplicate (frame));
751  cpl_frameset_join(preproc_frame, cpl_frameset_duplicate(used_frameset)); // add the calibration frameset
752 
753  applist = gravi_propertylist_get_qc (gravi_data_get_propertylist
754  (preproc_data, GRAVI_PRIMARY_HDR_NAME_EXT));
755 
756  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, PREPROC);
757 
758  gravi_data_save(preproc_data, frameset, preproc_name, parlist,
759  preproc_frame, frame, "gravi_dual", applist);
760 
761  CPLCHECK_CLEAN("Cannot save the PREPROC product");
762 
763  FREE (cpl_free,preproc_name);
764  FREE (cpl_propertylist_delete,applist);
765  FREE (cpl_frameset_delete,preproc_frame);
766  }
767 
768  /* Compute the flux and visibilities for each telescope and
769  * per acquisition with the P2VM applied to preproc_data */
770  cpl_msg_info (cpl_func, "Application of the p2vm to file %s", FILESHORT(filename));
771  p2vm_reduce = gravi_p2vm_reduce (preproc_data, p2vm_map, "gravi_dual", parlist);
772  FREE (gravi_data_delete,preproc_data);
773 
774  CPLCHECK_CLEAN("Cannot apply p2vm to the preproc data");
775 
776  /*
777  * Compute the QC0 about tau0 from piezo signals
778  */
779  gravi_compute_tau0 (p2vm_reduce);
780 
781  /* Visibility and flux are averaged and the followings
782  * are saved in Visibility data in tables VIS, VIS2 and T3 */
783  cpl_msg_info (cpl_func, "Computing the average visibilities of file %s", FILESHORT(filename));
784  gravi_vis_reduce (oi_vis, p2vm_reduce, nb_vis, parlist, "gravi_dual");
785 
786  CPLCHECK_CLEAN("Cannot average the visibilities");
787 
788  /* Compute the QC parameters of the TF
789  * TODO: to be done only if a calibration star */
790  cpl_msg_info( cpl_func, "Compute the QC TF parameters" );
791  cpl_errorstate prestate = cpl_errorstate_get();
792 
793  gravi_data * oi_tf = gravi_compute_tf( oi_vis );
794 
795  cpl_propertylist_copy_property_regexp( gravi_data_get_propertylist (oi_vis, GRAVI_PRIMARY_HDR_NAME_EXT),
796  gravi_data_get_propertylist (oi_tf, GRAVI_PRIMARY_HDR_NAME_EXT),
797  ".*QC TF.*", 0);
798  FREE (gravi_data_delete,oi_tf);
799 
800  /* If an error is catch when computing the QC parameters, dump this error but continue */
801  if( cpl_error_get_code() ) {
802  cpl_msg_warning(cpl_func, "Cannot compute the QC TF parameters for this observation");
803  cpl_errorstate_dump (prestate, 0, NULL);
804  cpl_error_reset();
805  }
806 
807  /* Save the p2vmreduced file */
808  if (cpl_parameter_get_bool(p_)){
809  cpl_msg_info( cpl_func, "Save the p2vmreduced file" );
810 
811  // char * p2vm_name = cpl_sprintf("dual_reduced.%s", FILESHORT(filename));
812  char * p2vm_name = gravi_data_product_name (filename, "p2vmreduced");
813 
814  cpl_frameset * p2vmreduced_frame = cpl_frameset_new ();
815  cpl_frameset_insert (p2vmreduced_frame, cpl_frame_duplicate (frame));
816  cpl_frameset_join(p2vmreduced_frame, cpl_frameset_duplicate(used_frameset));
817  applist = gravi_propertylist_get_qc (gravi_data_get_propertylist
818  (p2vm_reduce, GRAVI_PRIMARY_HDR_NAME_EXT));
819 
820  if ((strcmp(frame_tag, DUAL_CAL_SKY_RAW) == 0) ||
821  (strcmp(frame_tag, DUAL_SCI_SKY_RAW) == 0))
822  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, SKYREDUCED);
823  else
824  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, P2VMREDUCED);
825 
826  gravi_data_save(p2vm_reduce, frameset, p2vm_name, parlist,
827  p2vmreduced_frame, frame, "gravi_dual", applist);
828 
829  CPLCHECK_CLEAN("Cannot save the P2VMREDUCED product");
830 
831  FREE (cpl_free,p2vm_name);
832  FREE (cpl_propertylist_delete,applist);
833  FREE (cpl_frameset_delete,p2vmreduced_frame);
834  }
835 
836  nb_vis++;
837  cpl_msg_info(cpl_func,"Free the p2vmreduced");
838  FREE (gravi_data_delete,p2vm_reduce);
839  }
840  /* End loop on the input files to reduce */
841 
842  /* Eventually performn the flat with VIS_FLAT if provided */
843  if (gravi_param_get_bool(parlist, "gravi.vis_flat")) {
844  if (vis_flat)
845  {
846  cpl_msg_info (cpl_func, "Perform the visibility flat with the provided VIS_FLAT file");
847  gravi_flat_vis (oi_vis, vis_flat);
848  CPLCHECK_CLEAN("Cannot flat the OI_VIS");
849  }
850  else
851  {
852  cpl_msg_warning (cpl_func,"vis-flat option requested but no VIS_FLAT in the frameset");
853  }
854  }
855 
856  /* Save the output data file */
857  frame = cpl_frameset_get_position (recipe_frameset, 0);
858  frame_tag = cpl_frame_get_tag (frame);
859  primary_hdr = gravi_data_get_propertylist(oi_vis,
860  GRAVI_PRIMARY_HDR_NAME_EXT);
861  applist = gravi_propertylist_get_qc (primary_hdr);
862 
863  cpl_propertylist_append_int(applist, GRAVI_NIGHT_OBS,
864  cpl_propertylist_get_int (primary_hdr, GRAVI_NIGHT_OBS));
865  cpl_propertylist_append_string(applist, "DATE-OBS",
866  cpl_propertylist_get_string (primary_hdr, "DATE-OBS"));
867 
868 
869  if ((strcmp(frame_tag, GRAVI_DUAL_CALIB) == 0)) {
870  vis_name = gravi_data_product_name (filename, "dualcal");
871  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_DUAL_CALIB);
872  }
873  else if ((strcmp(frame_tag, GRAVI_DUAL_SCIENCE) == 0)) {
874  vis_name = gravi_data_product_name (filename, "dualsci");
875  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_DUAL_SCIENCE);
876  }
877  else {
878  vis_name = gravi_data_product_name (filename, "dualsci");
879  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_DUAL_SCIENCE);
880  }
881 
882  cpl_frameset_join (recipe_frameset, cpl_frameset_duplicate (used_frameset));
883 
884  gravi_data_save(oi_vis, frameset, vis_name, parlist,
885  recipe_frameset, frame, "gravi_dual", applist);
886 
887  CPLCHECK_CLEAN("Cannot save the VIS product");
888 
889 
890  /* Terminate the function */
891  goto cleanup;
892 
893 cleanup:
894  /* Deallocation of all variables */
895  cpl_msg_info(cpl_func,"Memory cleanup");
896 
897  FREE (gravi_data_delete,dark_map);
898  FREE (gravi_data_delete,data);
899  FREE (cpl_free,shutter);
900  FREE (cpl_free,preproc_name);
901  FREE (gravi_data_delete,preproc_data);
902  FREE (gravi_data_delete,profile_map);
903  FREE (gravi_data_delete,darkOrSky_map);
904  FREE (gravi_data_delete,wave_map);
905  FREE (gravi_data_delete,badpix_map);
906  FREE (gravi_data_delete,p2vm_map);
907  FREE (gravi_data_delete,vis_flat);
908  FREE (cpl_frameset_delete,darkcalib_frameset);
909  FREE (cpl_frameset_delete,wave_frameset);
910  FREE (cpl_frameset_delete,dark_frameset);
911  FREE (cpl_frameset_delete,gain_frameset);
912  FREE (cpl_frameset_delete,badpix_frameset);
913  FREE (gravi_data_delete,p2vm_reduce);
914  FREE (cpl_free,calib_datas);
915  FREE (cpl_free,vis_name);
916  FREE (cpl_propertylist_delete,applist);
917  FREE (gravi_data_delete,oi_vis);
918  FREE (cpl_frameset_delete,p2vm_frameset);
919  FREE (cpl_frameset_delete,visflat_frameset);
920  FREE (cpl_frameset_delete,recipe_frameset);
921  FREE (cpl_frameset_delete,used_frameset);
922 
923  cpl_msg_info(cpl_func,"Exit function");
924  return (int)cpl_error_get_code();
925 }
926 
927