GRAVI Pipeline Reference Manual  0.7.12
gravi_single.c
1 /* $Id: gravi_single.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 //#include <omp.h>
47 
48 
49 
50 /*-----------------------------------------------------------------------------
51  Private function prototypes
52  -----------------------------------------------------------------------------*/
53 
54 static int gravi_single_create(cpl_plugin *);
55 static int gravi_single_exec(cpl_plugin *);
56 static int gravi_single_destroy(cpl_plugin *);
57 static int gravi_single(cpl_frameset *, const cpl_parameterlist *);
58 
59 /*-----------------------------------------------------------------------------
60  Static variables
61  -----------------------------------------------------------------------------*/
62 
63 static char gravi_single_description[] =
64 "This recipe is associated to the observing template in single field mode.\n"
65 "Its aim is to reduce the raw interferometric data acquired on calibrator or \n"
66 "science targets and to compute the intrumental visibility.\n"
67 "The output OIFITS file contains the values of the complex visibilities, \n"
68 "squared visibilities and the closure phases.\n"
69 "Description DO category\n"
70 "Required input :\n"
71 " Master dark " DARK "\n"
72 " Master flat " FLAT "\n"
73 " Bad pixel map " BAD "\n"
74 " Master wave " WAVE "\n"
75 " Master p2vm " P2VM "\n"
76 " Raw observation file " GRAVI_SINGLE_CALIB " or "
77 GRAVI_SINGLE_SCIENCE "\n"
78 "Optional input\n"
79 " Sky file (one or more) " SINGLE_SCI_SKY_RAW " or "
80 SINGLE_CAL_SKY_RAW "\n"
81 
82 "Output :\n"
83 " Master sky (one per input sky) " SKY "\n"
84 " Raw Visibilities " VIS_SINGLE_CALIB " or " VIS_SINGLE_SCIENCE "\n"
85 "Optional output\n"
86 " Preprocessed files (--preproc-file) " PREPROC "\n"
87 " Non averaged visibilities (--p2vmreduced-file) " P2VMREDUCED "\n"
88 "\n";
89 
90 /*-----------------------------------------------------------------------------
91  Function code
92  -----------------------------------------------------------------------------*/
93 
94 /*----------------------------------------------------------------------------*/
104 /*----------------------------------------------------------------------------*/
105 int cpl_plugin_get_info(cpl_pluginlist * list)
106 {
107  cpl_recipe * recipe = cpl_calloc(1, sizeof *recipe );
108  cpl_plugin * plugin = &recipe->interface;
109 
110  if (cpl_plugin_init(plugin,
111  CPL_PLUGIN_API,
112  GRAVI_BINARY_VERSION,
113  CPL_PLUGIN_TYPE_RECIPE,
114  "gravi_single",
115  "This recipe is used to compute the squared visibility",
116  gravi_single_description,
117  "Firstname Lastname",
118  PACKAGE_BUGREPORT,
119  gravi_get_license(),
120  gravi_single_create,
121  gravi_single_exec,
122  gravi_single_destroy)) {
123  cpl_msg_error(cpl_func, "Plugin initialization failed");
124  (void)cpl_error_set_where(cpl_func);
125  return 1;
126  }
127 
128  if (cpl_pluginlist_append(list, plugin)) {
129  cpl_msg_error(cpl_func, "Error adding plugin to list");
130  (void)cpl_error_set_where(cpl_func);
131  return 1;
132  }
133 
134  return 0;
135 }
136 
137 /*----------------------------------------------------------------------------*/
145 /*----------------------------------------------------------------------------*/
146 static int gravi_single_create(cpl_plugin * plugin)
147 {
148  cpl_recipe * recipe;
149  cpl_parameter * p;
150 
151  /* Do not create the recipe if an error code is already set */
152  if (cpl_error_get_code() != CPL_ERROR_NONE) {
153  cpl_msg_error(cpl_func, "%s():%d: An error is already set: %s",
154  cpl_func, __LINE__, cpl_error_get_where());
155  return (int)cpl_error_get_code();
156  }
157 
158  if (plugin == NULL) {
159  cpl_msg_error(cpl_func, "Null plugin");
160  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
161  }
162 
163  /* Verify plugin type */
164  if (cpl_plugin_get_type(plugin) != CPL_PLUGIN_TYPE_RECIPE) {
165  cpl_msg_error(cpl_func, "Plugin is not a recipe");
166  cpl_ensure_code(0, (int)CPL_ERROR_TYPE_MISMATCH);
167  }
168 
169  /* Get the recipe */
170  recipe = (cpl_recipe *)plugin;
171 
172  /* Create the parameters list in the cpl_recipe object */
173  recipe->parameters = cpl_parameterlist_new();
174  if (recipe->parameters == NULL) {
175  cpl_msg_error(cpl_func, "Parameter list allocation failed");
176  cpl_ensure_code(0, (int)CPL_ERROR_ILLEGAL_OUTPUT);
177  }
178 
179  /* Fill the parameters list */
180 
181  /* --Save the preproc files */
182  p = cpl_parameter_new_value("gravi.preproc_param.preproc_file",
183  CPL_TYPE_BOOL, "Save the preprocessed file", "gravi.preproc_param", FALSE);
184  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "preproc-file");
185  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
186  cpl_parameterlist_append(recipe->parameters, p);
187 
188  /* --Save the p2vm reduced files */
189  p = cpl_parameter_new_value("gravi.preproc_param.p2vmreduced_file",
190  CPL_TYPE_BOOL, "Save the file with non averaged visibilities", "gravi.preproc_param", FALSE);
191  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "p2vmreduced-file");
192  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
193  cpl_parameterlist_append(recipe->parameters, p);
194 
195  /* -- Visibility correction */
196  p = cpl_parameter_new_value("gravi.vis_reduce.vis_correction",
197  CPL_TYPE_BOOL, "Correction of visibility losses", "gravi.vis_reduce", FALSE);
198  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "vis-correction");
199  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
200  cpl_parameterlist_append(recipe->parameters, p);
201 
202  /* --SNR threshold for fringe DET in FT */
203  p = cpl_parameter_new_value("gravi.snr_threshold_ft",
204  CPL_TYPE_DOUBLE, "SNR threshold for fringe detection in FT (>0)", "gravi.vis_reduce", 3.0);
205  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "snr-threshold-ft");
206  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
207  cpl_parameterlist_append(recipe->parameters, p);
208 
209  /* --STATE threshold for fringe DET in FT */
210  p = cpl_parameter_new_value("gravi.state_threshold_ft",
211  CPL_TYPE_DOUBLE, "STATE threshold for fringe detection in FT (>=0)", "gravi.vis_reduce", 1.0);
212  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "state-threshold-ft");
213  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
214  cpl_parameterlist_append(recipe->parameters, p);
215 
216  /* --Minimum detection ratio to accept SC frame */
217  p = cpl_parameter_new_value("gravi.fringedet_threshold_sc",
218  CPL_TYPE_DOUBLE, "Fringe-detection ratio threshold to accept SC frame (0..1)", "gravi.vis_reduce", 0.8);
219  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "fringedet-threshold-sc");
220  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
221  cpl_parameterlist_append(recipe->parameters, p);
222 
223  /* --vFactor threshold to accept SC frame */
224  p = cpl_parameter_new_value("gravi.vfactor_threshold_sc",
225  CPL_TYPE_DOUBLE, "vFactor threshold to accept SC frame (0..1)", "gravi.vis_reduce", 0.1);
226  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "vfactor-threshold-sc");
227  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
228  cpl_parameterlist_append(recipe->parameters, p);
229 
230  /* --Debias the V2 of SC */
231  p = cpl_parameter_new_value("gravi.debias_sc",
232  CPL_TYPE_BOOL, "Subtract the V2 bias from SC", "gravi", TRUE);
233  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "debias-sc");
234  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
235  cpl_parameterlist_append(recipe->parameters, p);
236 
237  /* --Debias the V2 of FT */
238  p = cpl_parameter_new_value("gravi.debias_ft",
239  CPL_TYPE_BOOL, "Subtract the V2 bias from FT", "gravi", TRUE);
240  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "debias-ft");
241  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
242  cpl_parameterlist_append(recipe->parameters, p);
243 
244  /* --Use the metrology */
245  p = cpl_parameter_new_value("gravi.use_met",
246  CPL_TYPE_BOOL, "Use the metrology to rephase the SC data in real-time", "gravi", FALSE);
247  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "use-met");
248  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
249  cpl_parameterlist_append(recipe->parameters, p);
250 
251  /* --Use the vis_flat */
252  p = cpl_parameter_new_value("gravi.vis_flat",
253  CPL_TYPE_BOOL, "Flat the OIVIS by the VIS_FLAT file", "gravi", TRUE);
254  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "vis-flat");
255  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
256  cpl_parameterlist_append(recipe->parameters, p);
257 
258  return 0;
259 }
260 
261 /*----------------------------------------------------------------------------*/
267 /*----------------------------------------------------------------------------*/
268 static int gravi_single_exec(cpl_plugin * plugin)
269 {
270 
271  cpl_recipe * recipe;
272  int recipe_status;
273  cpl_errorstate initial_errorstate = cpl_errorstate_get();
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_single(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_single_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_single(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, *bad_primary_hdr=NULL;
374  cpl_table * detector_table=NULL;
375  cpl_frame * frame=NULL;
376  const cpl_parameter * p=NULL, * p_=NULL;
377  cpl_parameter * par=NULL;
378  char * out_dark=NULL, * preproc_name=NULL;
379  const char * frame_tag=NULL, * filename=NULL, * temp=NULL;
380  gravi_data * p2vm_map=NULL, * data=NULL, * darkOrSky_map=NULL, * wave_map=NULL, * dark_map=NULL,
381  * profile_map=NULL, * badpix_map=NULL, * vis_flat=NULL;
382  gravi_data ** calib_datas=NULL;
383  gravi_data * preproc_data=NULL;
384  gravi_data ** raw_data=NULL, * p2vm_reduce=NULL;
385  gravi_data * oi_vis=NULL;
386  darkOrSky_map = NULL, wave_map = NULL;
387  int nb_frame, i, j, nwave_ft, nwave_sc, n_region;
388  int * shutter=NULL;
389  int nb_vis;
390  int testDarkSky = 0;
391  char * vis_name=NULL;
392 
393  /* Message */
394  cpl_msg_set_time_on();
395  cpl_msg_set_component_on();
396  cpl_msg_info(cpl_func,"Start function -- cleanup");
397 
398 
399  /*
400  * Identify the RAW and CALIB frames in the input frameset
401  */
402 
403  cpl_ensure_code(gravi_dfs_set_groups(frameset) == CPL_ERROR_NONE,
404  cpl_error_get_code()) ;
405 
406  p2vm_frameset = gravi_frameset_extract_p2vm_file(frameset);
407 
408  /* - Extract a set of frame p2vm and set of frame dark wave */
409  recipe_frameset = gravi_frameset_extract_single_data(frameset);
410  darkcalib_frameset = gravi_frameset_extract_dark_file(frameset);
411  wave_frameset = gravi_frameset_extract_wave_file(frameset);
412  dark_frameset = gravi_frameset_extract_dark(frameset);
413  gain_frameset = gravi_frameset_extract_flat_file(frameset);
414  badpix_frameset = gravi_frameset_extract_badpix(frameset);
415  skyframeset = gravi_frameset_extract_single_sky_data (frameset);
416  visflat_frameset = gravi_frameset_extract_visflat(frameset);
417 
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 SINGLE 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_single", 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 = cpl_sprintf("master_sky.%s", FILESHORT(filename));
532  out_dark = gravi_data_product_name (filename, "sky");
533 
534  gravi_data_save(darkOrSky_map, frameset, out_dark, parlist,
535  skyframeset, frame, "gravi_single", applist);
536 
537  CPLCHECK_CLEAN("Could not save the darkOrSky_map");
538 
539  FREE (cpl_propertylist_delete,applist);
540  FREE (cpl_free,out_dark);
541  }
542  else if (testDarkSky) {
543  darkOrSky_map = dark_map;
544  dark_map = NULL;
545  testDarkSky = 0;
546  cpl_msg_info (cpl_func, "There is no sky in the frame set");
547  }
548  else {
549  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE, "No dark or sky frames in the frameset");
550  goto cleanup;
551  }
552 
553  FREE (cpl_frameset_delete,skyframeset);
554 
555  /*
556  * Identify the BAD in the input frameset
557  */
558 
559  if (!cpl_frameset_is_empty(badpix_frameset)){
560 
561  /* Load the BAD pixel map */
562  frame = cpl_frameset_get_position(badpix_frameset, 0);
563  filename = cpl_frame_get_filename(frame);
564 
565  cpl_msg_info (cpl_func, "File %s is a BAD pixel map", FILESHORT(filename));
566  badpix_map = gravi_data_load(filename);
567  }
568  else {
569  /* Compute the BAD pixel map from the DARK or SKY */
570  badpix_map = gravi_compute_badpix(darkOrSky_map, parlist);
571 
572  CPLCHECK_CLEAN("Could not compute the BAD pixel map");
573 
574  /* Save the BAD pixel map */
575  bad_primary_hdr = gravi_data_get_propertylist (badpix_map, GRAVI_PRIMARY_HDR_NAME_EXT);
576  applist = gravi_propertylist_get_qc (bad_primary_hdr);
577 
578  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, BAD);
579 
580  frame = cpl_frameset_get_position (dark_frameset, 0);
581 
582  gravi_data_save (badpix_map, frameset, "gravi_bad_map.fits", parlist,
583  dark_frameset, frame, "gravi_single", applist);
584 
585  CPLCHECK_CLEAN("Could not save the BAD pixel map");
586 
587  FREE (cpl_propertylist_delete, applist);
588  }
589 
590  /*
591  * Identify the FLAT in the input frameset
592  */
593 
594  frame = cpl_frameset_get_position(gain_frameset, 0);
595  filename = cpl_frame_get_filename(frame);
596 
597  cpl_msg_info (cpl_func, "File %s is a master FLAT already computed", FILESHORT(filename));
598  profile_map = gravi_data_load(filename);
599 
600  /*
601  * Identify the WAVE in the input frameset
602  */
603 
604  frame = cpl_frameset_get_position(wave_frameset, 0);
605  filename = cpl_frame_get_filename(frame);
606 
607  cpl_msg_info (cpl_func, "File %s is a master WAVE map already computed", FILESHORT(filename));
608  wave_map = gravi_data_load(filename);
609 
610  /* Check if exist in the input frameset */
611  if ((wave_map == NULL) || (profile_map == NULL)) {
612  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE, "No WAVE or FLAT frames in the frameset");
613  goto cleanup;
614  }
615 
616  /*
617  * Identify the P2VM in the input frameset
618  */
619 
620  frame = cpl_frameset_get_position(p2vm_frameset, 0);
621  filename = cpl_frame_get_filename(frame);
622 
623  cpl_msg_info (cpl_func, "File %s is a P2VM file already computed", FILESHORT(filename));
624  p2vm_map = gravi_data_load (filename);
625 
626  CPLCHECK_CLEAN("Error while loading the P2VM map");
627 
628  /*
629  * Identify the VIS_FLAT in the input frameset
630  */
631  if (!cpl_frameset_is_empty(visflat_frameset))
632  {
633  frame = cpl_frameset_get_position(visflat_frameset, 0);
634  filename = cpl_frame_get_filename(frame);
635  cpl_msg_info (cpl_func, "File %s is a VIS_FLAT already computed", FILESHORT(filename));
636  vis_flat = gravi_data_load(filename);
637  }
638 
639  /*
640  * Prepare all the input and calibration for the reduction
641  */
642 
643  /* Check if the data with the p2vm have the same spectral n */
644  detector_table = gravi_data_get_table (p2vm_map, GRAVI_IMAGING_DETECTOR_FT_EXT);
645  n_region = cpl_table_get_nrow(detector_table);
646  if (n_region > 24) {
647  wave_plist = gravi_data_get_oi_propertylist
648  (p2vm_map, GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_SC),
649  GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_SC,0,2));
650  nwave_sc = gravi_pfits_get_nwave (wave_plist);
651  wave_plist = gravi_data_get_oi_propertylist (p2vm_map,
652  GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_FT),
653  GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_FT,0,2));
654  nwave_ft = gravi_pfits_get_nwave (wave_plist);
655  }
656  else {
657  wave_plist = gravi_data_get_oi_propertylist
658  (p2vm_map, GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_SC), GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_SC,0,1));
659  nwave_sc = gravi_pfits_get_nwave (wave_plist);
660  wave_plist = gravi_data_get_oi_propertylist (p2vm_map,
661  GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_FT), GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_FT,0,1));
662  nwave_ft = gravi_pfits_get_nwave (wave_plist);
663  }
664 
665  primary_hdr = gravi_data_get_propertylist (wave_map, GRAVI_PRIMARY_HDR_NAME_EXT);
666 
667  CPLCHECK_CLEAN("Missed one the preproc parameters");
668 
669  /* Construction of the calibration data */
670  int nb_calfile;
671  if (testDarkSky)
672  nb_calfile = 5;
673  else
674  nb_calfile = 4;
675 
676  calib_datas = cpl_malloc(nb_calfile * sizeof(gravi_data *));
677 
678  calib_datas[0] = darkOrSky_map;
679  calib_datas[1] = wave_map;
680  calib_datas[2] = profile_map;
681  calib_datas[3] = badpix_map;
682  if (testDarkSky)
683  calib_datas[4] = dark_map;
684 
685  /* insert calibration frame into the used frameset */
686  used_frameset=cpl_frameset_new();
687  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(dark_frameset));
688  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(badpix_frameset));
689  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(wave_frameset));
690  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(gain_frameset));
691  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(p2vm_frameset));
692 
693  p = cpl_parameterlist_find_const(parlist, "gravi.preproc_param.preproc_file");
694 
695  p_ = cpl_parameterlist_find_const(parlist, "gravi.preproc_param.p2vmreduced_file");
696 
697  /* Construction of the preproc data and extract the spectra using
698  * the results of the spectral calibration */
699 
700  nb_vis = 0;
701  oi_vis = gravi_data_new(0);
702  primary_hdr = gravi_data_get_propertylist (p2vm_map,
703  GRAVI_PRIMARY_HDR_NAME_EXT);
704 
705  const char * resol_p2vm = gravi_pfits_get_resolution (primary_hdr);
706  const char * resol_file;
707 
708 
709  /*
710  * Loop on input RAW frames to be reduced
711  */
712 
713  int thread;
714  for (thread = 0; thread < nb_frame; thread++){
715 
716  frame = cpl_frameset_get_position (recipe_frameset, thread);
717  frame_tag = cpl_frame_get_tag (frame);
718  filename = cpl_frame_get_filename (frame);
719 
720  data = gravi_data_load(filename);
721  primary_hdr = gravi_data_get_propertylist (data,
722  GRAVI_PRIMARY_HDR_NAME_EXT);
723 
724  /* Check the shutter for extracting the files can be used to compute
725  * the visibilities */
726  shutter = gravi_shutters_check(primary_hdr);
727  resol_file = gravi_pfits_get_resolution (primary_hdr);
728 
729  /* Check shutters */
730  if ((shutter [0] != 1) || (shutter [1] != 1) ||
731  (shutter [2] != 1) || (shutter [3] != 1)) {
732  cpl_msg_warning( cpl_func, "RAW data shutters: %d-%d-%d-%d",
733  shutter[0], shutter[1], shutter[2], shutter[3] );
734  }
735 
736  cpl_msg_info (cpl_func, "Preprocessing of the single"
737  " field file %s", FILESHORT(filename));
738 
739  preproc_data = gravi_preproc(data, calib_datas, nb_calfile, p2vm_map, parlist);
740 
741  FREE (gravi_data_delete,data);
742  FREE (cpl_free,shutter);
743 
744  CPLCHECK_CLEAN("Cannot preproc the data");
745 
746  /* Option save the proproc file */
747  if (cpl_parameter_get_bool(p)){
748  // preproc_name = cpl_sprintf("gravi_preproc.%s", FILESHORT(filename));
749  preproc_name = gravi_data_product_name (filename, "preproc");
750  cpl_frameset * preproc_frame = cpl_frameset_new ();
751  cpl_frameset_insert (preproc_frame, cpl_frame_duplicate (frame));
752  cpl_frameset_join(preproc_frame, cpl_frameset_duplicate(used_frameset)); // add the calibration frameset
753 
754  applist = gravi_propertylist_get_qc (gravi_data_get_propertylist
755  (preproc_data, GRAVI_PRIMARY_HDR_NAME_EXT));
756 
757  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, PREPROC);
758 
759  gravi_data_save(preproc_data, frameset, preproc_name, parlist,
760  preproc_frame, frame, "gravi_single", applist);
761 
762  CPLCHECK_CLEAN("Cannot save the PREPROC product");
763 
764  FREE (cpl_free,preproc_name);
765  FREE (cpl_propertylist_delete,applist);
766  FREE (cpl_frameset_delete,preproc_frame);
767  }
768 
769  /* Compute the flux and visibilities for each telescope and
770  * per acquisition with the P2VM applied to preproc_data */
771  cpl_msg_info (cpl_func, "Application of the p2vm to file %s", FILESHORT(filename));
772  p2vm_reduce = gravi_p2vm_reduce (preproc_data, p2vm_map, "gravi_single", parlist);
773  FREE (gravi_data_delete,preproc_data);
774 
775  CPLCHECK_CLEAN("Cannot apply p2vm to the preproc data");
776 
777  /*
778  * Compute the QC0 about tau0 from piezo signals
779  */
780  gravi_compute_tau0 (p2vm_reduce);
781 
782  /* Visibility and flux are averaged and the followings
783  * are saved in Visibility data in tables VIS, VIS2 and T3 */
784  cpl_msg_info (cpl_func, "Computing the average visibilities of file %s", FILESHORT(filename));
785  gravi_vis_reduce (oi_vis, p2vm_reduce, nb_vis, parlist, "gravi_single");
786 
787  CPLCHECK_CLEAN("Cannot average the visibilities");
788 
789  /* Compute the QC parameters of the TF
790  * TODO: to be done only if a calibration star */
791  cpl_msg_info( cpl_func, "Compute the QC TF parameters" );
792  cpl_errorstate prestate = cpl_errorstate_get();
793 
794  gravi_data * oi_tf = gravi_compute_tf( oi_vis );
795 
796  cpl_propertylist_copy_property_regexp( gravi_data_get_propertylist (oi_vis, GRAVI_PRIMARY_HDR_NAME_EXT),
797  gravi_data_get_propertylist (oi_tf, GRAVI_PRIMARY_HDR_NAME_EXT),
798  ".*QC TF.*", 0);
799  FREE (gravi_data_delete,oi_tf);
800 
801  /* If an error is catch when computing the QC parameters, dump this error but continue */
802  if( cpl_error_get_code() ) {
803  cpl_msg_warning(cpl_func, "Cannot compute the QC TF parameters for this observation");
804  cpl_errorstate_dump (prestate, 0, NULL);
805  cpl_error_reset();
806  }
807 
808  /* Save the p2vmreduced file */
809  if (cpl_parameter_get_bool(p_)){
810  cpl_msg_info( cpl_func, "Save the p2vmreduced file" );
811 
812  // char * p2vm_name = cpl_sprintf("single_reduced.%s", FILESHORT(filename));
813  char * p2vm_name = gravi_data_product_name (filename, "p2vmreduced");
814 
815  cpl_frameset * p2vmreduced_frame = cpl_frameset_new ();
816  cpl_frameset_insert (p2vmreduced_frame, cpl_frame_duplicate (frame));
817  cpl_frameset_join(p2vmreduced_frame, cpl_frameset_duplicate(used_frameset));
818  applist = gravi_propertylist_get_qc (gravi_data_get_propertylist
819  (p2vm_reduce, GRAVI_PRIMARY_HDR_NAME_EXT));
820 
821  if ((strcmp(frame_tag, SINGLE_CAL_SKY_RAW) == 0) ||
822  (strcmp(frame_tag, SINGLE_SCI_SKY_RAW) == 0))
823  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, SKYREDUCED);
824  else
825  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, P2VMREDUCED);
826 
827  gravi_data_save(p2vm_reduce, frameset, p2vm_name, parlist,
828  p2vmreduced_frame, frame, "gravi_single", applist);
829 
830  CPLCHECK_CLEAN("Cannot save the P2VMREDUCED product");
831 
832  FREE (cpl_free,p2vm_name);
833  FREE (cpl_propertylist_delete,applist);
834  FREE (cpl_frameset_delete,p2vmreduced_frame);
835  }
836 
837  nb_vis++;
838  cpl_msg_info(cpl_func,"Free the p2vmreduced");
839  FREE (gravi_data_delete,p2vm_reduce);
840  }
841  /* End loop on the input files to reduce */
842 
843 
844  /* Eventually performn the flat with VIS_FLAT if provided */
845  if (gravi_param_get_bool(parlist, "gravi.vis_flat")) {
846  if (vis_flat)
847  {
848  cpl_msg_info (cpl_func, "Perform the visibility flat with the provided VIS_FLAT file");
849  gravi_flat_vis (oi_vis, vis_flat);
850  CPLCHECK_CLEAN("Cannot flat the OI_VIS");
851  }
852  else
853  {
854  cpl_msg_warning (cpl_func,"vis-flat option requested but no VIS_FLAT in the frameset");
855  }
856  }
857 
858 
859 
860  /* Save the output data file */
861  frame = cpl_frameset_get_position (recipe_frameset, 0);
862  frame_tag = cpl_frame_get_tag (frame);
863  primary_hdr = gravi_data_get_propertylist(oi_vis,
864  GRAVI_PRIMARY_HDR_NAME_EXT);
865  applist = gravi_propertylist_get_qc (primary_hdr);
866 
867  cpl_propertylist_append_int(applist, GRAVI_NIGHT_OBS,
868  cpl_propertylist_get_int (primary_hdr, GRAVI_NIGHT_OBS));
869  cpl_propertylist_append_string(applist, "DATE-OBS",
870  cpl_propertylist_get_string (primary_hdr, "DATE-OBS"));
871 
872 
873  if ((strcmp(frame_tag, GRAVI_SINGLE_CALIB) == 0)) {
874  vis_name = gravi_data_product_name (filename, "singlecal");
875  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_SINGLE_CALIB);
876  }
877  else if ((strcmp(frame_tag, GRAVI_SINGLE_SCIENCE) == 0)) {
878  vis_name = gravi_data_product_name (filename, "singlesci");
879  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_SINGLE_SCIENCE);
880  }
881  else {
882  vis_name = gravi_data_product_name (filename, "singlesci");
883  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_SINGLE_SCIENCE);
884  }
885 
886  cpl_frameset_join (recipe_frameset, cpl_frameset_duplicate (used_frameset));
887 
888  gravi_data_save(oi_vis, frameset, vis_name, parlist,
889  recipe_frameset, frame, "gravi_single", applist);
890 
891  CPLCHECK_CLEAN("Cannot save the VIS product");
892 
893 
894  /* Terminate the function */
895  goto cleanup;
896 
897 cleanup:
898  /* Deallocation of all variables */
899  cpl_msg_info(cpl_func,"Memory cleanup");
900 
901  FREE (gravi_data_delete,dark_map);
902  FREE (gravi_data_delete,data);
903  FREE (cpl_free,shutter);
904  FREE (cpl_free,preproc_name);
905  FREE (gravi_data_delete,preproc_data);
906  FREE (gravi_data_delete,profile_map);
907  FREE (gravi_data_delete,darkOrSky_map);
908  FREE (gravi_data_delete,wave_map);
909  FREE (gravi_data_delete,badpix_map);
910  FREE (gravi_data_delete,p2vm_map);
911  FREE (gravi_data_delete,vis_flat);
912  FREE (cpl_frameset_delete,darkcalib_frameset);
913  FREE (cpl_frameset_delete,wave_frameset);
914  FREE (cpl_frameset_delete,dark_frameset);
915  FREE (cpl_frameset_delete,gain_frameset);
916  FREE (cpl_frameset_delete,badpix_frameset);
917  FREE (gravi_data_delete,p2vm_reduce);
918  FREE (cpl_free,calib_datas);
919  FREE (cpl_free,vis_name);
920  FREE (cpl_propertylist_delete,applist);
921  FREE (gravi_data_delete,oi_vis);
922  FREE (cpl_frameset_delete,p2vm_frameset);
923  FREE (cpl_frameset_delete,visflat_frameset);
924  FREE (cpl_frameset_delete,recipe_frameset);
925  FREE (cpl_frameset_delete,used_frameset);
926 
927  cpl_msg_info(cpl_func,"Exit function");
928  return (int)cpl_error_get_code();
929 }
930 
931