GRAVI Pipeline Reference Manual  0.7.11
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 visibilities looses", "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  /* --Minimum locking ratio to accept SC frame */
210  p = cpl_parameter_new_value("gravi.fringedet_threshold_sc",
211  CPL_TYPE_DOUBLE, "Fringe-detection ratio threshold to accept SC frame (0..1)", "gravi.vis_reduce", 0.8);
212  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "fringedet-threshold-sc");
213  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
214  cpl_parameterlist_append(recipe->parameters, p);
215 
216  /* --Debias the V2 of SC */
217  p = cpl_parameter_new_value("gravi.debias_sc",
218  CPL_TYPE_BOOL, "Subtrast the V2 bias from SC", "gravi", TRUE);
219  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "debias-sc");
220  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
221  cpl_parameterlist_append(recipe->parameters, p);
222 
223  /* --Debias the V2 of FT */
224  p = cpl_parameter_new_value("gravi.debias_ft",
225  CPL_TYPE_BOOL, "Subtrast the V2 bias from FT", "gravi", TRUE);
226  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "debias-ft");
227  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
228  cpl_parameterlist_append(recipe->parameters, p);
229 
230  /* --Use the metrology */
231  p = cpl_parameter_new_value("gravi.use_met",
232  CPL_TYPE_BOOL, "Use the metrology to rephase the SC data in real-time", "gravi", FALSE);
233  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "use-met");
234  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
235  cpl_parameterlist_append(recipe->parameters, p);
236 
237  /* --Use the vis_flat */
238  p = cpl_parameter_new_value("gravi.vis_flat",
239  CPL_TYPE_BOOL, "Flat the OIVIS by the VIS_FLAT file", "gravi", TRUE);
240  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "vis-flat");
241  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
242  cpl_parameterlist_append(recipe->parameters, p);
243 
244  return 0;
245 }
246 
247 /*----------------------------------------------------------------------------*/
253 /*----------------------------------------------------------------------------*/
254 static int gravi_single_exec(cpl_plugin * plugin)
255 {
256 
257  cpl_recipe * recipe;
258  int recipe_status;
259  cpl_errorstate initial_errorstate = cpl_errorstate_get();
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_single(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_single_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_single(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, *bad_primary_hdr=NULL;
360  cpl_table * detector_table=NULL;
361  cpl_frame * frame=NULL;
362  const cpl_parameter * p=NULL, * p_=NULL;
363  cpl_parameter * par=NULL;
364  char * out_dark=NULL, * preproc_name=NULL;
365  const char * frame_tag=NULL, * filename=NULL, * temp=NULL;
366  gravi_data * p2vm_map=NULL, * data=NULL, * darkOrSky_map=NULL, * wave_map=NULL, * dark_map=NULL,
367  * profile_map=NULL, * badpix_map=NULL, * vis_flat=NULL;
368  gravi_data ** calib_datas=NULL;
369  gravi_data * preproc_data=NULL;
370  gravi_data ** raw_data=NULL, * p2vm_reduce=NULL;
371  gravi_data * oi_vis=NULL;
372  darkOrSky_map = NULL, wave_map = NULL;
373  int nb_frame, i, j, nwave_ft, nwave_sc, n_region;
374  int * shutter=NULL;
375  int nb_vis;
376  int testDarkSky = 0;
377  char * vis_name=NULL;
378 
379  /* Message */
380  cpl_msg_set_time_on();
381  cpl_msg_set_component_on();
382  cpl_msg_info(cpl_func,"Start function -- cleanup");
383 
384 
385  /*
386  * Identify the RAW and CALIB frames in the input frameset
387  */
388 
389  cpl_ensure_code(gravi_dfs_set_groups(frameset) == CPL_ERROR_NONE,
390  cpl_error_get_code()) ;
391 
392  p2vm_frameset = gravi_frameset_extract_p2vm_file(frameset);
393 
394  /* - Extract a set of frame p2vm and set of frame dark wave */
395  recipe_frameset = gravi_frameset_extract_single_data(frameset);
396  darkcalib_frameset = gravi_frameset_extract_dark_file(frameset);
397  wave_frameset = gravi_frameset_extract_wave_file(frameset);
398  dark_frameset = gravi_frameset_extract_dark(frameset);
399  gain_frameset = gravi_frameset_extract_flat_file(frameset);
400  badpix_frameset = gravi_frameset_extract_badpix(frameset);
401  skyframeset = gravi_frameset_extract_single_sky_data (frameset);
402  visflat_frameset = gravi_frameset_extract_visflat(frameset);
403 
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 SINGLE 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_single", 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 = cpl_sprintf("master_sky.%s", FILESHORT(filename));
518  out_dark = gravi_data_product_name (filename, "sky");
519 
520  gravi_data_save(darkOrSky_map, frameset, out_dark, parlist,
521  skyframeset, frame, "gravi_single", applist);
522 
523  CPLCHECK_CLEAN("Could not save the darkOrSky_map");
524 
525  FREE (cpl_propertylist_delete,applist);
526  FREE (cpl_free,out_dark);
527  }
528  else if (testDarkSky) {
529  darkOrSky_map = dark_map;
530  dark_map = NULL;
531  testDarkSky = 0;
532  cpl_msg_info (cpl_func, "There is no sky in the frame set");
533  }
534  else {
535  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE, "No dark or sky frames in the frameset");
536  goto cleanup;
537  }
538 
539  FREE (cpl_frameset_delete,skyframeset);
540 
541  /*
542  * Identify the BAD in the input frameset
543  */
544 
545  if (!cpl_frameset_is_empty(badpix_frameset)){
546 
547  /* Load the BAD pixel map */
548  frame = cpl_frameset_get_position(badpix_frameset, 0);
549  filename = cpl_frame_get_filename(frame);
550 
551  cpl_msg_info (cpl_func, "File %s is a BAD pixel map", FILESHORT(filename));
552  badpix_map = gravi_data_load(filename);
553  }
554  else {
555  /* Compute the BAD pixel map from the DARK or SKY */
556  badpix_map = gravi_compute_badpix(darkOrSky_map, parlist);
557 
558  CPLCHECK_CLEAN("Could not compute the BAD pixel map");
559 
560  /* Save the BAD pixel map */
561  bad_primary_hdr = gravi_data_get_propertylist (badpix_map, GRAVI_PRIMARY_HDR_NAME_EXT);
562  applist = gravi_propertylist_get_qc (bad_primary_hdr);
563 
564  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, BAD);
565 
566  frame = cpl_frameset_get_position (dark_frameset, 0);
567 
568  gravi_data_save (badpix_map, frameset, "gravi_bad_map.fits", parlist,
569  dark_frameset, frame, "gravi_single", applist);
570 
571  CPLCHECK_CLEAN("Could not save the BAD pixel map");
572 
573  FREE (cpl_propertylist_delete, applist);
574  }
575 
576  /*
577  * Identify the FLAT in the input frameset
578  */
579 
580  frame = cpl_frameset_get_position(gain_frameset, 0);
581  filename = cpl_frame_get_filename(frame);
582 
583  cpl_msg_info (cpl_func, "File %s is a master FLAT already computed", FILESHORT(filename));
584  profile_map = gravi_data_load(filename);
585 
586  /*
587  * Identify the WAVE in the input frameset
588  */
589 
590  frame = cpl_frameset_get_position(wave_frameset, 0);
591  filename = cpl_frame_get_filename(frame);
592 
593  cpl_msg_info (cpl_func, "File %s is a master WAVE map already computed", FILESHORT(filename));
594  wave_map = gravi_data_load(filename);
595 
596  /* Check if exist in the input frameset */
597  if ((wave_map == NULL) || (profile_map == NULL)) {
598  cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE, "No WAVE or FLAT frames in the frameset");
599  goto cleanup;
600  }
601 
602  /*
603  * Identify the P2VM in the input frameset
604  */
605 
606  frame = cpl_frameset_get_position(p2vm_frameset, 0);
607  filename = cpl_frame_get_filename(frame);
608 
609  cpl_msg_info (cpl_func, "File %s is a P2VM file already computed", FILESHORT(filename));
610  p2vm_map = gravi_data_load (filename);
611 
612  CPLCHECK_CLEAN("Error while loading the P2VM map");
613 
614  /*
615  * Identify the VIS_FLAT in the input frameset
616  */
617  if (!cpl_frameset_is_empty(visflat_frameset))
618  {
619  frame = cpl_frameset_get_position(visflat_frameset, 0);
620  filename = cpl_frame_get_filename(frame);
621  cpl_msg_info (cpl_func, "File %s is a VIS_FLAT already computed", FILESHORT(filename));
622  vis_flat = gravi_data_load(filename);
623  }
624 
625  /*
626  * Prepare all the input and calibration for the reduction
627  */
628 
629  /* Check if the data with the p2vm have the same spectral n */
630  detector_table = gravi_data_get_table (p2vm_map, GRAVI_IMAGING_DETECTOR_FT_EXT);
631  n_region = cpl_table_get_nrow(detector_table);
632  if (n_region > 24) {
633  wave_plist = gravi_data_get_oi_propertylist
634  (p2vm_map, GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_SC),
635  GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_SC,0,2));
636  nwave_sc = gravi_pfits_get_nwave (wave_plist);
637  wave_plist = gravi_data_get_oi_propertylist (p2vm_map,
638  GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_FT),
639  GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_FT,0,2));
640  nwave_ft = gravi_pfits_get_nwave (wave_plist);
641  }
642  else {
643  wave_plist = gravi_data_get_oi_propertylist
644  (p2vm_map, GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_SC), GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_SC,0,1));
645  nwave_sc = gravi_pfits_get_nwave (wave_plist);
646  wave_plist = gravi_data_get_oi_propertylist (p2vm_map,
647  GRAVI_OI_WAVE_SWITCH(GRAVI_TYPE_FT), GRAVI_SPECTRO_SWITCH(GRAVI_TYPE_FT,0,1));
648  nwave_ft = gravi_pfits_get_nwave (wave_plist);
649  }
650 
651  primary_hdr = gravi_data_get_propertylist (wave_map, GRAVI_PRIMARY_HDR_NAME_EXT);
652 
653  CPLCHECK_CLEAN("Missed one the preproc parameters");
654 
655  /* Construction of the calibration data */
656  int nb_calfile;
657  if (testDarkSky)
658  nb_calfile = 5;
659  else
660  nb_calfile = 4;
661 
662  calib_datas = cpl_malloc(nb_calfile * sizeof(gravi_data *));
663 
664  calib_datas[0] = darkOrSky_map;
665  calib_datas[1] = wave_map;
666  calib_datas[2] = profile_map;
667  calib_datas[3] = badpix_map;
668  if (testDarkSky)
669  calib_datas[4] = dark_map;
670 
671  /* insert calibration frame into the used frameset */
672  used_frameset=cpl_frameset_new();
673  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(dark_frameset));
674  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(badpix_frameset));
675  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(wave_frameset));
676  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(gain_frameset));
677  cpl_frameset_join(used_frameset, cpl_frameset_duplicate(p2vm_frameset));
678 
679  p = cpl_parameterlist_find_const(parlist, "gravi.preproc_param.preproc_file");
680 
681  p_ = cpl_parameterlist_find_const(parlist, "gravi.preproc_param.p2vmreduced_file");
682 
683  /* Construction of the preproc data and extract the spectra using
684  * the results of the spectral calibration */
685 
686  nb_vis = 0;
687  oi_vis = gravi_data_new(0);
688  primary_hdr = gravi_data_get_propertylist (p2vm_map,
689  GRAVI_PRIMARY_HDR_NAME_EXT);
690 
691  const char * resol_p2vm = gravi_pfits_get_resolution (primary_hdr);
692  const char * resol_file;
693 
694 
695  /*
696  * Loop on input RAW frames to be reduced
697  */
698 
699  int thread;
700  for (thread = 0; thread < nb_frame; thread++){
701 
702  frame = cpl_frameset_get_position (recipe_frameset, thread);
703  frame_tag = cpl_frame_get_tag (frame);
704  filename = cpl_frame_get_filename (frame);
705 
706  data = gravi_data_load(filename);
707  primary_hdr = gravi_data_get_propertylist (data,
708  GRAVI_PRIMARY_HDR_NAME_EXT);
709 
710  /* Check the shutter for extracting the files can be used to compute
711  * the visibilities */
712  shutter = gravi_shutters_check(primary_hdr);
713  resol_file = gravi_pfits_get_resolution (primary_hdr);
714 
715  /* Check shutters */
716  if ((shutter [0] != 1) || (shutter [1] != 1) ||
717  (shutter [2] != 1) || (shutter [3] != 1)) {
718  cpl_msg_warning( cpl_func, "RAW data shutters: %d-%d-%d-%d",
719  shutter[0], shutter[1], shutter[2], shutter[3] );
720  }
721 
722  cpl_msg_info (cpl_func, "Preprocessing of the single"
723  " field file %s", FILESHORT(filename));
724 
725  preproc_data = gravi_preproc(data, calib_datas, nb_calfile, p2vm_map, parlist);
726 
727  FREE (gravi_data_delete,data);
728  FREE (cpl_free,shutter);
729 
730  CPLCHECK_CLEAN("Cannot preproc the data");
731 
732  /* Option save the proproc file */
733  if (cpl_parameter_get_bool(p)){
734  // preproc_name = cpl_sprintf("gravi_preproc.%s", FILESHORT(filename));
735  preproc_name = gravi_data_product_name (filename, "preproc");
736  cpl_frameset * preproc_frame = cpl_frameset_new ();
737  cpl_frameset_insert (preproc_frame, cpl_frame_duplicate (frame));
738  cpl_frameset_join(preproc_frame, cpl_frameset_duplicate(used_frameset)); // add the calibration frameset
739 
740  applist = gravi_propertylist_get_qc (gravi_data_get_propertylist
741  (preproc_data, GRAVI_PRIMARY_HDR_NAME_EXT));
742 
743  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, PREPROC);
744 
745  gravi_data_save(preproc_data, frameset, preproc_name, parlist,
746  preproc_frame, frame, "gravi_single", applist);
747 
748  CPLCHECK_CLEAN("Cannot save the PREPROC product");
749 
750  FREE (cpl_free,preproc_name);
751  FREE (cpl_propertylist_delete,applist);
752  FREE (cpl_frameset_delete,preproc_frame);
753  }
754 
755  /* Compute the flux and visibilities for each telescope and
756  * per acquisition with the P2VM applied to preproc_data */
757  cpl_msg_info (cpl_func, "Application of the p2vm to file %s", FILESHORT(filename));
758  p2vm_reduce = gravi_p2vm_reduce (preproc_data, p2vm_map, "gravi_single", parlist);
759  FREE (gravi_data_delete,preproc_data);
760 
761  CPLCHECK_CLEAN("Cannot apply p2vm to the preproc data");
762 
763  /*
764  * Compute the QC0 about tau0 from piezo signals
765  */
766  gravi_compute_tau0 (p2vm_reduce);
767 
768  /* Visibility and flux are averaged and the followings
769  * are saved in Visibility data in tables VIS, VIS2 and T3 */
770  cpl_msg_info (cpl_func, "Computing the average visibilities of file %s", FILESHORT(filename));
771  gravi_vis_reduce (oi_vis, p2vm_reduce, nb_vis, parlist, "gravi_single");
772 
773  CPLCHECK_CLEAN("Cannot average the visibilities");
774 
775  /* Compute the QC parameters of the TF
776  * TODO: to be done only if a calibration star */
777  cpl_msg_info( cpl_func, "Compute the QC TF parameters" );
778  cpl_errorstate prestate = cpl_errorstate_get();
779 
780  gravi_data * oi_tf = gravi_compute_tf( oi_vis );
781 
782  cpl_propertylist_copy_property_regexp( gravi_data_get_propertylist (oi_vis, GRAVI_PRIMARY_HDR_NAME_EXT),
783  gravi_data_get_propertylist (oi_tf, GRAVI_PRIMARY_HDR_NAME_EXT),
784  ".*QC TF.*", 0);
785  FREE (gravi_data_delete,oi_tf);
786 
787  /* If an error is catch when computing the QC parameters, dump this error but continue */
788  if( cpl_error_get_code() ) {
789  cpl_msg_warning(cpl_func, "Cannot compute the QC TF parameters for this observation");
790  cpl_errorstate_dump (prestate, 0, NULL);
791  cpl_error_reset();
792  }
793 
794  /* Save the p2vmreduced file */
795  if (cpl_parameter_get_bool(p_)){
796  cpl_msg_info( cpl_func, "Save the p2vmreduced file" );
797 
798  // char * p2vm_name = cpl_sprintf("single_reduced.%s", FILESHORT(filename));
799  char * p2vm_name = gravi_data_product_name (filename, "p2vmreduced");
800 
801  cpl_frameset * p2vmreduced_frame = cpl_frameset_new ();
802  cpl_frameset_insert (p2vmreduced_frame, cpl_frame_duplicate (frame));
803  cpl_frameset_join(p2vmreduced_frame, cpl_frameset_duplicate(used_frameset));
804  applist = gravi_propertylist_get_qc (gravi_data_get_propertylist
805  (p2vm_reduce, GRAVI_PRIMARY_HDR_NAME_EXT));
806 
807  if ((strcmp(frame_tag, SINGLE_CAL_SKY_RAW) == 0) ||
808  (strcmp(frame_tag, SINGLE_SCI_SKY_RAW) == 0))
809  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, SKYREDUCED);
810  else
811  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, P2VMREDUCED);
812 
813  gravi_data_save(p2vm_reduce, frameset, p2vm_name, parlist,
814  p2vmreduced_frame, frame, "gravi_single", applist);
815 
816  CPLCHECK_CLEAN("Cannot save the P2VMREDUCED product");
817 
818  FREE (cpl_free,p2vm_name);
819  FREE (cpl_propertylist_delete,applist);
820  FREE (cpl_frameset_delete,p2vmreduced_frame);
821  }
822 
823  nb_vis++;
824  cpl_msg_info(cpl_func,"Free the p2vmreduced");
825  FREE (gravi_data_delete,p2vm_reduce);
826  }
827  /* End loop on the input files to reduce */
828 
829 
830  /* Eventually performn the flat with VIS_FLAT if provided */
831  if (gravi_param_get_bool(parlist, "gravi.vis_flat")) {
832  if (vis_flat)
833  {
834  cpl_msg_info (cpl_func, "Perform the visibility flat with the provided VIS_FLAT file");
835  gravi_flat_vis (oi_vis, vis_flat);
836  CPLCHECK_CLEAN("Cannot flat the OI_VIS");
837  }
838  else
839  {
840  cpl_msg_warning (cpl_func,"vis-flat option requested but no VIS_FLAT in the frameset");
841  }
842  }
843 
844 
845 
846  /* Save the output data file */
847  frame = cpl_frameset_get_position (recipe_frameset, 0);
848  frame_tag = cpl_frame_get_tag (frame);
849  primary_hdr = gravi_data_get_propertylist(oi_vis,
850  GRAVI_PRIMARY_HDR_NAME_EXT);
851  applist = gravi_propertylist_get_qc (primary_hdr);
852 
853  cpl_propertylist_append_int(applist, GRAVI_NIGHT_OBS,
854  cpl_propertylist_get_int (primary_hdr, GRAVI_NIGHT_OBS));
855  cpl_propertylist_append_string(applist, "DATE-OBS",
856  cpl_propertylist_get_string (primary_hdr, "DATE-OBS"));
857 
858 
859  if ((strcmp(frame_tag, GRAVI_SINGLE_CALIB) == 0)) {
860  vis_name = gravi_data_product_name (filename, "singlecal");
861  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_SINGLE_CALIB);
862  }
863  else if ((strcmp(frame_tag, GRAVI_SINGLE_SCIENCE) == 0)) {
864  vis_name = gravi_data_product_name (filename, "singlesci");
865  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_SINGLE_SCIENCE);
866  }
867  else {
868  vis_name = gravi_data_product_name (filename, "singlesci");
869  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_SINGLE_SCIENCE);
870  }
871 
872  cpl_frameset_join (recipe_frameset, cpl_frameset_duplicate (used_frameset));
873 
874  gravi_data_save(oi_vis, frameset, vis_name, parlist,
875  recipe_frameset, frame, "gravi_single", applist);
876 
877  CPLCHECK_CLEAN("Cannot save the VIS product");
878 
879 
880  /* Terminate the function */
881  goto cleanup;
882 
883 cleanup:
884  /* Deallocation of all variables */
885  cpl_msg_info(cpl_func,"Memory cleanup");
886 
887  FREE (gravi_data_delete,dark_map);
888  FREE (gravi_data_delete,data);
889  FREE (cpl_free,shutter);
890  FREE (cpl_free,preproc_name);
891  FREE (gravi_data_delete,preproc_data);
892  FREE (gravi_data_delete,profile_map);
893  FREE (gravi_data_delete,darkOrSky_map);
894  FREE (gravi_data_delete,wave_map);
895  FREE (gravi_data_delete,badpix_map);
896  FREE (gravi_data_delete,p2vm_map);
897  FREE (gravi_data_delete,vis_flat);
898  FREE (cpl_frameset_delete,darkcalib_frameset);
899  FREE (cpl_frameset_delete,wave_frameset);
900  FREE (cpl_frameset_delete,dark_frameset);
901  FREE (cpl_frameset_delete,gain_frameset);
902  FREE (cpl_frameset_delete,badpix_frameset);
903  FREE (gravi_data_delete,p2vm_reduce);
904  FREE (cpl_free,calib_datas);
905  FREE (cpl_free,vis_name);
906  FREE (cpl_propertylist_delete,applist);
907  FREE (gravi_data_delete,oi_vis);
908  FREE (cpl_frameset_delete,p2vm_frameset);
909  FREE (cpl_frameset_delete,visflat_frameset);
910  FREE (cpl_frameset_delete,recipe_frameset);
911  FREE (cpl_frameset_delete,used_frameset);
912 
913  cpl_msg_info(cpl_func,"Exit function");
914  return (int)cpl_error_get_code();
915 }
916 
917