GRAVI Pipeline Reference Manual  0.6.3
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 
47 
48 
49 /*-----------------------------------------------------------------------------
50  Private function prototypes
51  -----------------------------------------------------------------------------*/
52 
53 static int gravi_single_create(cpl_plugin *);
54 static int gravi_single_exec(cpl_plugin *);
55 static int gravi_single_destroy(cpl_plugin *);
56 static int gravi_single(cpl_frameset *, const cpl_parameterlist *);
57 
58 /*-----------------------------------------------------------------------------
59  Static variables
60  -----------------------------------------------------------------------------*/
61 
62 static char gravi_single_description[] =
63 "This recipe is associated to the template GRAVI_single_calib.\n"
64 "Its aim is to reduce the raw interferometric data acquired on calibrator\n"
65 "target and to compute the intrumental visibility. It is used in single mode.\n"
66 "The given output FITS file contain a oi fits file, \n"
67 "the values of the visibility complex, squared visibility and the cloture phase.\n"
68 "their associated tags, e.g."
69 "GRAVI-gravi_single-raw-file.fits " GRAVI_SINGLE_CALIB " and "
70 GRAVI_SINGLE_SCIENCE "\n"
71 "and any optional files, e.g.\n"
72 "GRAVI-gravi_single-file.fits " VIS_CALIB " and " VIS_SCIENCE "\n"
73 "Additionally, it should describe functionality of the expected output."
74 "\n";
75 
76 /*-----------------------------------------------------------------------------
77  Function code
78  -----------------------------------------------------------------------------*/
79 
80 /*----------------------------------------------------------------------------*/
90 /*----------------------------------------------------------------------------*/
91 int cpl_plugin_get_info(cpl_pluginlist * list)
92 {
93  cpl_recipe * recipe = cpl_calloc(1, sizeof *recipe );
94  cpl_plugin * plugin = &recipe->interface;
95 
96  if (cpl_plugin_init(plugin,
97  CPL_PLUGIN_API,
98  GRAVI_BINARY_VERSION,
99  CPL_PLUGIN_TYPE_RECIPE,
100  "gravi_single",
101  "This recipe is used to compute the squared visibility",
102  gravi_single_description,
103  "Firstname Lastname",
104  PACKAGE_BUGREPORT,
105  gravi_get_license(),
106  gravi_single_create,
107  gravi_single_exec,
108  gravi_single_destroy)) {
109  cpl_msg_error(cpl_func, "Plugin initialization failed");
110  (void)cpl_error_set_where(cpl_func);
111  return 1;
112  }
113 
114  if (cpl_pluginlist_append(list, plugin)) {
115  cpl_msg_error(cpl_func, "Error adding plugin to list");
116  (void)cpl_error_set_where(cpl_func);
117  return 1;
118  }
119 
120  return 0;
121 }
122 
123 /*----------------------------------------------------------------------------*/
131 /*----------------------------------------------------------------------------*/
132 static int gravi_single_create(cpl_plugin * plugin)
133 {
134  cpl_recipe * recipe;
135  cpl_parameter * p;
136 
137  /* Do not create the recipe if an error code is already set */
138  if (cpl_error_get_code() != CPL_ERROR_NONE) {
139  cpl_msg_error(cpl_func, "%s():%d: An error is already set: %s",
140  cpl_func, __LINE__, cpl_error_get_where());
141  return (int)cpl_error_get_code();
142  }
143 
144  if (plugin == NULL) {
145  cpl_msg_error(cpl_func, "Null plugin");
146  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
147  }
148 
149  /* Verify plugin type */
150  if (cpl_plugin_get_type(plugin) != CPL_PLUGIN_TYPE_RECIPE) {
151  cpl_msg_error(cpl_func, "Plugin is not a recipe");
152  cpl_ensure_code(0, (int)CPL_ERROR_TYPE_MISMATCH);
153  }
154 
155  /* Get the recipe */
156  recipe = (cpl_recipe *)plugin;
157 
158  /* Create the parameters list in the cpl_recipe object */
159  recipe->parameters = cpl_parameterlist_new();
160  if (recipe->parameters == NULL) {
161  cpl_msg_error(cpl_func, "Parameter list allocation failed");
162  cpl_ensure_code(0, (int)CPL_ERROR_ILLEGAL_OUTPUT);
163  }
164 
165  /* Fill the parameters list */
166 
167  /* --the sigma readout */
168  p = cpl_parameter_new_value("gravi.gravi_all_p2vm.preproc_param.readout",
169  CPL_TYPE_DOUBLE, "the readout error sigma", "gravi.gravi_all_p2vm", 1);
170  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "sigma_readout");
171  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
172  cpl_parameterlist_append(recipe->parameters, p);
173 
174  /* --Save the preproc files */
175  p = cpl_parameter_new_value("gravi.gravi_all_p2vm.preproc_param.preproc_file",
176  CPL_TYPE_BOOL, "Save the preproc file", "gravi.gravi_all_p2vm", FALSE);
177  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "preprocfile");
178  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
179  cpl_parameterlist_append(recipe->parameters, p);
180 
181  /* --Save the preproc files */
182  p = cpl_parameter_new_value("gravi.gravi_all_p2vm.vis_reduce.vis_correction",
183  CPL_TYPE_BOOL, "Parameter of corection", "gravi.gravi_single", FALSE);
184  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "vis_correction");
185  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
186  cpl_parameterlist_append(recipe->parameters, p);
187 
188  /* --the sigma readout SC */
189  p = cpl_parameter_new_value("gravi.gravi_all_p2vm.preproc_param.readoutSC",
190  CPL_TYPE_DOUBLE, "the readout error sigma SC", "gravi.gravi_all_p2vm", 1.0);
191  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "sigma_readout_SC");
192  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
193  cpl_parameterlist_append(recipe->parameters, p);
194 
195  /* --the sigma readout FT */
196  p = cpl_parameter_new_value("gravi.gravi_all_p2vm.preproc_param.readoutFT",
197  CPL_TYPE_DOUBLE, "the readout error sigma FT", "gravi.gravi_all_p2vm", 1.0);
198  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "sigma_readout_FT");
199  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
200  cpl_parameterlist_append(recipe->parameters, p);
201 
202  return 0;
203 }
204 
205 /*----------------------------------------------------------------------------*/
211 /*----------------------------------------------------------------------------*/
212 static int gravi_single_exec(cpl_plugin * plugin)
213 {
214 
215  cpl_recipe * recipe;
216  int recipe_status;
217  cpl_errorstate initial_errorstate = cpl_errorstate_get();
218 
219  /* Return immediately if an error code is already set */
220  if (cpl_error_get_code() != CPL_ERROR_NONE) {
221  cpl_msg_error(cpl_func, "%s():%d: An error is already set: %s",
222  cpl_func, __LINE__, cpl_error_get_where());
223  return (int)cpl_error_get_code();
224  }
225 
226  if (plugin == NULL) {
227  cpl_msg_error(cpl_func, "Null plugin");
228  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
229  }
230 
231  /* Verify plugin type */
232  if (cpl_plugin_get_type(plugin) != CPL_PLUGIN_TYPE_RECIPE) {
233  cpl_msg_error(cpl_func, "Plugin is not a recipe");
234  cpl_ensure_code(0, (int)CPL_ERROR_TYPE_MISMATCH);
235  }
236 
237  /* Get the recipe */
238  recipe = (cpl_recipe *)plugin;
239 
240  /* Verify parameter and frame lists */
241  if (recipe->parameters == NULL) {
242  cpl_msg_error(cpl_func, "Recipe invoked with NULL parameter list");
243  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
244  }
245  if (recipe->frames == NULL) {
246  cpl_msg_error(cpl_func, "Recipe invoked with NULL frame set");
247  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
248  }
249 
250  /* Invoke the recipe */
251  recipe_status = gravi_single(recipe->frames, recipe->parameters);
252 
253  /* Ensure DFS-compliance of the products */
254 
255  if (cpl_dfs_update_product_header(recipe->frames)) {
256  if (!recipe_status){
257  recipe_status = (int)cpl_error_get_code();
258  }
259  }
260 
261  if (!cpl_errorstate_is_equal(initial_errorstate)) {
262  /* Dump the error history since recipe execution start.
263  At this point the recipe cannot recover from the error */
264  cpl_errorstate_dump(initial_errorstate, CPL_FALSE, NULL);
265  }
266 
267  return recipe_status;
268 }
269 
270 /*----------------------------------------------------------------------------*/
276 /*----------------------------------------------------------------------------*/
277 static int gravi_single_destroy(cpl_plugin * plugin)
278 {
279  cpl_recipe * recipe;
280 
281  if (plugin == NULL) {
282  cpl_msg_error(cpl_func, "Null plugin");
283  cpl_ensure_code(0, (int)CPL_ERROR_NULL_INPUT);
284  }
285 
286  /* Verify plugin type */
287  if (cpl_plugin_get_type(plugin) != CPL_PLUGIN_TYPE_RECIPE) {
288  cpl_msg_error(cpl_func, "Plugin is not a recipe");
289  cpl_ensure_code(0, (int)CPL_ERROR_TYPE_MISMATCH);
290  }
291 
292  /* Get the recipe */
293  recipe = (cpl_recipe *)plugin;
294 
295  cpl_parameterlist_delete(recipe->parameters);
296 
297  return 0;
298 }
299 
300 
301 /*----------------------------------------------------------------------------*/
309 /*----------------------------------------------------------------------------*/
310 static int gravi_single(cpl_frameset * frameset,
311  const cpl_parameterlist * parlist)
312 {
313  cpl_frameset * recipe_frameset, * wave_frameset, * dark_frameset,
314  * darkcalib_frameset, * preproc_frameset, * skyframeset,
315  * gain_frameset, * p2vm_frameset, * badpix_frameset;
316  cpl_propertylist * applist, * primary_hdr, * primary_hdr_dark, * wave_plist
317  , *bad_primary_hdr;
318 
319  cpl_parameterlist * preproc_parlist;
320  cpl_frame * frame;
321  const cpl_parameter * p;
322  cpl_parameter * par;
323  char * out_dark, * preproc_name;
324  const char * frame_tag, * filename/*, * insname*/;
325  gravi_data * p2vm_map, * data, * dark_map, * wave_map,
326  * profile_map, * badpix_map;
327  gravi_data ** calib_datas;
328  gravi_data * preproc_data;
329  gravi_data ** raw_data, * p2vm_reduce;
330  gravi_data * oi_vis;
331  dark_map = NULL, wave_map = NULL;
332  int nb_frame, i, j, nwave_ft, nwave_sc;
333  int * shutter;
334  int nb_vis;
335  char * vis_name;
336  double darkmean, darkrms, wavemin, wavemax;
337 
338  /* Identify the RAW and CALIB frames in the input frameset */
339 
340  cpl_ensure_code(gravi_dfs_set_groups(frameset) == CPL_ERROR_NONE,
341  cpl_error_get_code()) ;
342 
343  preproc_frameset = gravi_frameset_extract_preproc(frameset);
344  p2vm_frameset = gravi_frameset_extract_p2vm_file(frameset);
345  if (cpl_frameset_is_empty(preproc_frameset)) {
346 
347  /* - Extract a set of frame p2vm and set of frame dark wave */
348 // cpl_frameset * frameset = cpl_frameset_duplicate(_frameset);
349  recipe_frameset = gravi_frameset_extract_single_data(frameset);
350 // p2vm_frameset = gravi_frameset_extract_p2vm_file(frameset);
351  darkcalib_frameset = gravi_frameset_extract_dark_file(frameset);
352  wave_frameset = gravi_frameset_extract_wave_file(frameset);
353  dark_frameset = gravi_frameset_extract_dark(frameset);
354  gain_frameset = gravi_frameset_extract_flat_file(frameset);
355  badpix_frameset = gravi_frameset_extract_badpix(frameset);
356  skyframeset = gravi_frameset_extract_sky_data (frameset);
357  if (cpl_frameset_is_empty(recipe_frameset)) {
358  /* To use this recipe the frameset must contain at least
359  * one p2vm frame. */
360  cpl_frameset_delete(darkcalib_frameset);
361  cpl_frameset_delete(wave_frameset);
362  cpl_frameset_delete(dark_frameset);
363  cpl_frameset_delete(gain_frameset);
364  cpl_frameset_delete(p2vm_frameset);
365  cpl_frameset_delete(recipe_frameset);
366 
367  return (int)cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE,
368  "No single frame on the frameset") ;
369  }
370 
371  if ((cpl_frameset_is_empty(p2vm_frameset)) ||
372  (cpl_frameset_is_empty(wave_frameset)) ||
373  (cpl_frameset_is_empty(gain_frameset))) {
374  /* To use this recipe the frameset must contain the p2vm, wave and
375  * gain calibration file. */
376  cpl_frameset_delete(darkcalib_frameset);
377  cpl_frameset_delete(recipe_frameset);
378  cpl_frameset_delete(wave_frameset);
379  cpl_frameset_delete(dark_frameset);
380  cpl_frameset_delete(gain_frameset);
381  cpl_frameset_delete(p2vm_frameset);
382  return (int)cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE,
383  "No p2vm or gain or wave file on the frameset") ;
384  }
385 
386  /* Get the number of the p2vm frame contained in the frameset */
387  nb_frame = cpl_frameset_get_size(recipe_frameset);
388 
389  /* Identify and extract the dark and wave fits files */
390 
391  /* Check and get the dark frame */
392 
393  if (cpl_frameset_is_empty(skyframeset)) {
394 
395  if (!cpl_frameset_is_empty(dark_frameset)) {
396  frame = cpl_frameset_get_position(dark_frameset, 0);
397 
398  filename = cpl_frame_get_filename(frame);
399  data = gravi_data_load(filename);
400  primary_hdr = gravi_data_get_propertylist(data,
401  GRAVI_PRIMARY_HDR_NAME_EXT);
402  shutter = gravi_shutters_check(primary_hdr);
403  if (!((shutter[0] == 0) && (shutter[1] == 0) && (shutter[2] == 0)
404  && (shutter[3] == 0))){
405  cpl_frameset_delete(darkcalib_frameset);
406  cpl_frameset_delete(dark_frameset);
407  cpl_frameset_delete(wave_frameset);
408  cpl_frameset_delete(gain_frameset);
409  cpl_frameset_delete(p2vm_frameset);
410  cpl_frameset_delete(recipe_frameset);
411  gravi_data_delete(data);
412  return (int)cpl_error_set_message(cpl_func,
413  CPL_ERROR_ILLEGAL_OUTPUT, "Shutter problem");
414  }
415  /* The dark file must contains all the shutter close */
416 
417  cpl_msg_info (NULL, "This file %s is a dark file",
418  filename);
419  dark_map = gravi_compute_dark(data);
420  if (cpl_error_get_code()) {
421  cpl_frameset_delete(darkcalib_frameset);
422  cpl_frameset_delete(dark_frameset);
423  cpl_frameset_delete(wave_frameset);
424  cpl_frameset_delete(gain_frameset);
425  gravi_data_delete(dark_map);
426  cpl_frameset_delete(p2vm_frameset);
427  cpl_frameset_delete(recipe_frameset);
428  cpl_propertylist_delete(applist);
429  gravi_data_delete(data);
430  return (int)cpl_error_set_message(cpl_func,
431  CPL_ERROR_ILLEGAL_OUTPUT,
432  "Error while computing the dark_map");
433  }
434 
435  /* Save the dark map */
436 
437  primary_hdr_dark = gravi_data_get_propertylist(dark_map,
438  GRAVI_PRIMARY_HDR_NAME_EXT);
439  applist = gravi_propertylist_get_qc (primary_hdr_dark);
440 // dit = gravi_pfits_get_dit(primary_hdr_dark);
441  //insname = gravi_pfits_get_insname(primary_hdr_dark);
442  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG,
443  DARK);
444 
445 // cpl_propertylist_append_double(applist, GRAVI_DET_DIT, dit);
446  //cpl_propertylist_append_string(applist, "INSNAME", insname);
447  cpl_propertylist_append_double(applist, QC_MEANDARK, darkmean);
448  cpl_propertylist_append_double(applist, QC_DARKRMS, darkrms);
449  out_dark = "gravi_final_dark.fits";
450 
451  if (gravi_data_save(dark_map, frameset, out_dark, parlist,
452  dark_frameset, frame, "gravi_single", applist)
453  != CPL_ERROR_NONE){
454  cpl_frameset_delete(darkcalib_frameset);
455  cpl_frameset_delete(dark_frameset);
456  cpl_frameset_delete(wave_frameset);
457  cpl_frameset_delete(gain_frameset);
458  gravi_data_delete(dark_map);
459  cpl_frameset_delete(p2vm_frameset);
460  cpl_frameset_delete(recipe_frameset);
461  cpl_propertylist_delete(applist);
462  gravi_data_delete(data);
463  cpl_free(shutter);
464 
465  return (int) cpl_error_set_message(cpl_func,
466  CPL_ERROR_ILLEGAL_OUTPUT, "Could not save the dark_map"
467  " on the output file");
468  }
469  cpl_propertylist_delete(applist);
470 
471 
472  gravi_data_delete(data);
473  cpl_free(shutter);
474 
475 
476  }
477  else if (!cpl_frameset_is_empty(darkcalib_frameset)) {
478  frame = cpl_frameset_get_position(darkcalib_frameset, 0);
479 
480  frame_tag = cpl_frame_get_tag(frame) ;
481  filename = cpl_frame_get_filename(frame);
482  cpl_msg_info (NULL, "This file %s is a dark file already "
483  "computed", filename);
484  dark_map = gravi_data_load(filename);
485  cpl_frameset_insert (dark_frameset, cpl_frame_duplicate (frame));
486  }
487 
488  /* Check that the dark_map in the input frameset */
489  else {
490  cpl_frameset_delete(darkcalib_frameset);
491  cpl_frameset_delete(p2vm_frameset);
492  cpl_frameset_delete(dark_frameset);
493  cpl_frameset_delete(wave_frameset);
494  cpl_frameset_delete(gain_frameset);
495  gravi_data_delete(dark_map);
496  cpl_frameset_delete(recipe_frameset);
497 
498  return (int)cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE,
499  "No dark frames in the frameset");
500  }
501  }
502  else {
503  frame = cpl_frameset_get_position(skyframeset, 0);
504 
505  filename = cpl_frame_get_filename(frame);
506  data = gravi_data_load(filename);
507  primary_hdr = gravi_data_get_propertylist(data,
508  GRAVI_PRIMARY_HDR_NAME_EXT);
509  const char * dpr_TYPE = cpl_propertylist_get_string (primary_hdr, DPR_TYPE);
510  if (cpl_error_get_code()){
511  cpl_frameset_delete(darkcalib_frameset);
512  cpl_frameset_delete(dark_frameset);
513  cpl_frameset_delete(wave_frameset);
514  cpl_frameset_delete(gain_frameset);
515  cpl_frameset_delete(p2vm_frameset);
516  cpl_frameset_delete(recipe_frameset);
517  cpl_frameset_delete(skyframeset);
518  gravi_data_delete(data);
519  return (int)cpl_error_set_message(cpl_func,
520  CPL_ERROR_ILLEGAL_OUTPUT, "The keyword DPR_CATG is missed");
521  }
522 
523  if (strstr (dpr_TYPE, "SKY") == NULL){
524  cpl_frameset_delete(darkcalib_frameset);
525  cpl_frameset_delete(dark_frameset);
526  cpl_frameset_delete(wave_frameset);
527  cpl_frameset_delete(gain_frameset);
528  cpl_frameset_delete(p2vm_frameset);
529  cpl_frameset_delete(recipe_frameset);
530  cpl_frameset_delete(skyframeset);
531  gravi_data_delete(data);
532  return (int)cpl_error_set_message(cpl_func,
533  CPL_ERROR_ILLEGAL_OUTPUT, "The file is not a sky file");
534  }
535  /* The dark file must contains all the shutter close */
536 
537  cpl_msg_info (NULL, "This file %s is a sky file",
538  filename);
539  dark_map = gravi_compute_dark(data);
540 
541  gravi_data_delete(data);
542  if (cpl_error_get_code()) {
543  cpl_frameset_delete(darkcalib_frameset);
544  cpl_frameset_delete(dark_frameset);
545  cpl_frameset_delete(wave_frameset);
546  cpl_frameset_delete(gain_frameset);
547  gravi_data_delete(dark_map);
548  cpl_frameset_delete(p2vm_frameset);
549  cpl_frameset_delete(recipe_frameset);
550  cpl_propertylist_delete(applist);
551 
552  return (int)cpl_error_set_message(cpl_func,
553  CPL_ERROR_ILLEGAL_OUTPUT,
554  "Error while computing the dark_map");
555  }
556  cpl_frameset_insert (dark_frameset, cpl_frame_duplicate (frame));
557  /* Save the dark map */
558 
559  primary_hdr_dark = gravi_data_get_propertylist(dark_map,
560  GRAVI_PRIMARY_HDR_NAME_EXT);
561  applist = gravi_propertylist_get_qc (primary_hdr_dark);
562 // dit = gravi_pfits_get_dit(primary_hdr_dark);
563  //insname = gravi_pfits_get_insname(primary_hdr_dark);
564  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG,
565  DARK);
566 
567 // cpl_propertylist_append_double(applist, GRAVI_DET_DIT, dit);
568  //cpl_propertylist_append_string(applist, "INSNAME", insname);
569  cpl_propertylist_append_double(applist, QC_MEANDARK, darkmean);
570  cpl_propertylist_append_double(applist, QC_DARKRMS, darkrms);
571  out_dark = "gravi_final_dark.fits";
572 
573  if (gravi_data_save(dark_map, frameset, out_dark, parlist,
574  skyframeset, frame, "gravi_single", applist)
575  != CPL_ERROR_NONE){
576  cpl_frameset_delete(darkcalib_frameset);
577  cpl_frameset_delete(dark_frameset);
578  cpl_frameset_delete(wave_frameset);
579  cpl_frameset_delete(gain_frameset);
580  gravi_data_delete(dark_map);
581  cpl_frameset_delete(p2vm_frameset);
582  cpl_frameset_delete(recipe_frameset);
583  cpl_propertylist_delete(applist);
584 
585  cpl_free(shutter);
586 
587  return (int) cpl_error_set_message(cpl_func,
588  CPL_ERROR_ILLEGAL_OUTPUT, "Could not save the dark_map"
589  " on the output file");
590  }
591  cpl_propertylist_delete(applist);
592  }
593 
594  cpl_frameset_delete(skyframeset);
595  /* Compute the bad pixel map */
596  if (!cpl_frameset_is_empty(badpix_frameset)){
597  frame = cpl_frameset_get_position(badpix_frameset, 0);
598  filename = cpl_frame_get_filename(frame);
599  cpl_msg_info (NULL, "This file %s is a bad pixel map", filename);
600 
601  badpix_map = gravi_data_load(filename);
602  }
603  else {
604 
605  badpix_map = gravi_compute_badpix(dark_map, parlist);
606 
607  if (cpl_error_get_code()) {
608  gravi_data_delete(badpix_map);
609  cpl_frameset_delete(dark_frameset);
610  cpl_frameset_delete(wave_frameset);
611  cpl_frameset_delete(gain_frameset);
612  cpl_frameset_delete(p2vm_frameset);
613  cpl_frameset_delete(darkcalib_frameset);
614  gravi_data_delete(dark_map);
615  return (int)cpl_error_set_message(cpl_func, CPL_ERROR_ILLEGAL_OUTPUT,
616  "Error while computing the bad pixel map");
617 
618  }
619  /* Save the pad pixel map */
620  bad_primary_hdr = gravi_data_get_propertylist (badpix_map,
621  GRAVI_PRIMARY_HDR_NAME_EXT);
622  applist = gravi_propertylist_get_qc (bad_primary_hdr);
623 
624  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, BAD);
625 
626  frame = cpl_frameset_get_position (dark_frameset, 0);
627 
628  if (gravi_data_save (badpix_map, frameset, "gravi_bad_map.fits", parlist,
629  dark_frameset, frame, "gravi_all_flat", applist)
630  != CPL_ERROR_NONE){
631 
632  gravi_data_delete(badpix_map);
633  cpl_frameset_delete(dark_frameset);
634  cpl_frameset_delete(wave_frameset);
635  cpl_frameset_delete(gain_frameset);
636  cpl_frameset_delete(p2vm_frameset);
637  cpl_frameset_delete(darkcalib_frameset);
638  cpl_frameset_delete(badpix_frameset);
639  gravi_data_delete(dark_map);
640  return (int) cpl_error_set_message(cpl_func,
641  CPL_ERROR_ILLEGAL_OUTPUT, "Could not save the bad pixel map");
642  }
643  cpl_propertylist_delete (applist);
644  }
645 
646  /* Identify the flat file */
647 
648  frame = cpl_frameset_get_position(gain_frameset, 0);
649  filename = cpl_frame_get_filename(frame);
650  cpl_msg_info (NULL, "This file %s is a flat file already "
651  "computed", filename);
652 
653  profile_map = gravi_data_load(filename);
654 
655  /* Check and get the wave frame */
656  frame = cpl_frameset_get_position(wave_frameset, 0);
657  filename = cpl_frame_get_filename(frame);
658  cpl_msg_info (NULL, "This file %s is a wave file already "
659  "computed", filename);
660  wave_map = gravi_data_load(filename);
661 
662  /* Check that wave_map exist in the input frameset */
663  if ((wave_map == NULL) || (profile_map == NULL)){
664  cpl_frameset_delete(darkcalib_frameset);
665  cpl_frameset_delete(p2vm_frameset);
666  cpl_frameset_delete(dark_frameset);
667  cpl_frameset_delete(wave_frameset);
668  cpl_frameset_delete(gain_frameset);
669  gravi_data_delete(wave_map);
670  gravi_data_delete(profile_map);
671  gravi_data_delete(badpix_map);
672  cpl_frameset_delete(recipe_frameset);
673 
674  return (int)cpl_error_set_message(cpl_func, CPL_ERROR_INVALID_TYPE,
675  "No wave frames in the frameset");
676  }
677 
678  /* Load the p2vm data */
679  frame = cpl_frameset_get_position(p2vm_frameset, 0);
680  filename = cpl_frame_get_filename(frame);
681  cpl_msg_info (NULL, "This file %s is a p2vm file already "
682  "computed", filename);
683 
684  p2vm_map = gravi_data_load (filename);
685 
686  if (cpl_error_get_code()) {
687  cpl_frameset_delete(darkcalib_frameset);
688  cpl_frameset_delete(p2vm_frameset);
689  cpl_frameset_delete(dark_frameset);
690  cpl_frameset_delete(wave_frameset);
691  cpl_frameset_delete(gain_frameset);
692  gravi_data_delete(profile_map);
693  gravi_data_delete(dark_map);
694  gravi_data_delete(wave_map);
695  gravi_data_delete(p2vm_map);
696  cpl_frameset_delete(recipe_frameset);
697 
698  return (int)cpl_error_set_message(cpl_func, CPL_ERROR_ILLEGAL_OUTPUT,
699  "Error while loading the p2vm_map");
700  }
701 
702  /* Check if the data with the p2vm have the same spectral n */
703  wave_plist = gravi_data_get_propertylist (p2vm_map,
704  GRAVI_IO_WAVELENGTH_SC_EXT);
705  nwave_sc = gravi_pfits_get_nwave (wave_plist);
706  wave_plist = gravi_data_get_propertylist (p2vm_map,
707  GRAVI_IO_WAVELENGTH_FT_EXT);
708  nwave_ft = gravi_pfits_get_nwave (wave_plist);
709  primary_hdr = gravi_data_get_propertylist (wave_map,
710  GRAVI_PRIMARY_HDR_NAME_EXT);//pour ft et sc
711 
712  if (cpl_error_get_code()) {
713  cpl_frameset_delete(darkcalib_frameset);
714  cpl_frameset_delete(p2vm_frameset);
715  cpl_frameset_delete(dark_frameset);
716  cpl_frameset_delete(wave_frameset);
717  cpl_frameset_delete(gain_frameset);
718  gravi_data_delete(profile_map);
719  gravi_data_delete(dark_map);
720  gravi_data_delete(wave_map);
721  gravi_data_delete(p2vm_map);
722  cpl_frameset_delete(recipe_frameset);
723 
724  return (int)cpl_error_set_message(cpl_func, CPL_ERROR_ILLEGAL_OUTPUT,
725  "Missed one the preproc parameters");
726  }
727 
728  wavemin = cpl_propertylist_get_double (primary_hdr, QC_MINWAVE_SC);
729  wavemax = cpl_propertylist_get_double (primary_hdr, QC_MAXWAVE_SC);
730 
731  /* Construction of the parameterlist to compute the preprocing */
732  cpl_size nb_par = cpl_parameterlist_get_size (parlist);
733 
734  preproc_parlist = cpl_parameterlist_new ();
735  p = cpl_parameterlist_get_first_const (parlist);
736  for (i = 0; i < nb_par; i ++){
737  cpl_type type = cpl_parameter_get_type (p);
738  switch (type) {
739  case CPL_TYPE_DOUBLE :
740  par = cpl_parameter_new_value(cpl_parameter_get_name (p), type,
741  cpl_parameter_get_help(p),
742  cpl_parameter_get_context(p),
743  cpl_parameter_get_double (p));
744  break;
745  case CPL_TYPE_BOOL :
746  par = cpl_parameter_new_value(cpl_parameter_get_name (p), type,
747  cpl_parameter_get_help(p),
748  cpl_parameter_get_context(p),
749  cpl_parameter_get_bool (p));
750  break;
751  case CPL_TYPE_INT :
752  par = cpl_parameter_new_value(cpl_parameter_get_name (p), type,
753  cpl_parameter_get_help(p),
754  cpl_parameter_get_context(p),
755  cpl_parameter_get_int (p));
756  break;
757  case CPL_TYPE_STRING :
758  par = cpl_parameter_new_value(cpl_parameter_get_name (p), type,
759  cpl_parameter_get_help(p),
760  cpl_parameter_get_context(p),
761  cpl_parameter_get_string (p));
762  break;
763  default :
764  cpl_msg_warning(cpl_func, "Type of the parameter is different %s",
765  cpl_type_get_name (type));
766  break;
767  }
768 
769 
770  cpl_parameter_set_alias(par, CPL_PARAMETER_MODE_CLI,
771  cpl_parameter_get_alias (p, CPL_PARAMETER_MODE_CLI));
772  cpl_parameter_disable(par, CPL_PARAMETER_MODE_ENV);
773  cpl_parameterlist_append(preproc_parlist, par);
774  p = cpl_parameterlist_get_next_const (parlist);
775  }
776  /* --the minimum wavelength */
777 
778  par = cpl_parameter_new_value("gravi.gravi_all_p2vm.preproc_param.min_wave",
779  CPL_TYPE_DOUBLE, "the minimum wavelength",
780  "gravi.gravi_all_p2vm", wavemin + 0.0001);
781  cpl_parameter_set_alias(par, CPL_PARAMETER_MODE_CLI, "lambdamin");
782  cpl_parameter_disable(par, CPL_PARAMETER_MODE_ENV);
783  cpl_parameterlist_append(preproc_parlist, par);
784 
785  /* --the maximum wavelength */
786  par = cpl_parameter_new_value("gravi.gravi_all_p2vm.preproc_param.max_wave",
787  CPL_TYPE_DOUBLE, "the maximum wavelength",
788  "gravi.gravi_all_p2vm", wavemax - 0.0001);
789  cpl_parameter_set_alias(par, CPL_PARAMETER_MODE_CLI, "lambdamax");
790  cpl_parameter_disable(par, CPL_PARAMETER_MODE_ENV);
791  cpl_parameterlist_append(preproc_parlist, par);
792 
793  /* --the number of the spectral element to compute */
794  par = cpl_parameter_new_value("gravi.gravi_all_p2vm.preproc_param.nspectrum_FT",
795  CPL_TYPE_INT, "the number of the spectral",
796  "gravi.gravi_all_dark", nwave_ft);
797  cpl_parameter_set_alias(par, CPL_PARAMETER_MODE_CLI, "nlambda_ft");
798  cpl_parameter_disable(par, CPL_PARAMETER_MODE_ENV);
799  cpl_parameterlist_append(preproc_parlist, par);
800 
801  /* --the number of the spectral element to compute */
802  par = cpl_parameter_new_value("gravi.gravi_all_p2vm.preproc_param.nspectrum_SC",
803  CPL_TYPE_INT, "the number of the spectral",
804  "gravi.gravi_all_dark", nwave_sc);
805  cpl_parameter_set_alias(par, CPL_PARAMETER_MODE_CLI, "nlambda_sc");
806  cpl_parameter_disable(par, CPL_PARAMETER_MODE_ENV);
807  cpl_parameterlist_append(preproc_parlist, par);
808 
809 
810 // printf("nb_frameset = %d nb_frameset_dark = %d nb_frameset_wave = %d nb_frameset_p2vm = %d nb_gain_frame = %d\n",
811 // nb_frame, cpl_frameset_get_size(dark_frameset), cpl_frameset_get_size(wave_frameset),
812 // cpl_frameset_get_size(p2vm_frameset), cpl_frameset_get_size(gain_frameset));
813 
814 
815  /* Construction of the calibration data */
816  calib_datas = cpl_malloc(4 * sizeof(gravi_data *));
817  calib_datas[0] = dark_map;
818  calib_datas[1] = wave_map;
819  calib_datas[2] = profile_map;
820  calib_datas[3] = badpix_map;
821 
822  p = cpl_parameterlist_find_const(parlist, "gravi.gravi_all_p2vm."
823  "preproc_param.preproc_file");
824 
825  /* Construction of the preproc data and extract the spectra using
826  * the results of the spectral calibration */
827 
828  nb_vis = 0;
829  oi_vis = gravi_data_new(0);
830  for (i = 0; i < nb_frame; i++){
831  frame = cpl_frameset_get_position (recipe_frameset, i);
832  filename = cpl_frame_get_filename (frame);
833 
834  data = gravi_data_load(filename);
835  primary_hdr = gravi_data_get_propertylist (data,
836  GRAVI_PRIMARY_HDR_NAME_EXT);
837  /* Check the shutter for extracting the files can be used to compute
838  * the visibilities */
839  shutter = gravi_shutters_check(primary_hdr);
840 
841  if ((shutter[0] == 1) && (shutter[1] == 1) &&
842  (shutter[2] == 1) && (shutter[3] == 1)){
843  cpl_msg_info (NULL, "Calibration of this file %s", filename);
844  preproc_data = gravi_preproc(data, calib_datas, 4, preproc_parlist);
845 
846  gravi_data_delete(data);
847  cpl_free(shutter);
848  cpl_frameset_insert (preproc_frameset, cpl_frame_duplicate (frame));
849  }
850  else {
851  gravi_data_delete(data);
852  cpl_free(shutter);
853  continue;
854  }
855  if (cpl_error_get_code()) {
856  cpl_frameset_delete(darkcalib_frameset);
857  cpl_frameset_delete(p2vm_frameset);
858  cpl_frameset_delete(dark_frameset);
859  cpl_frameset_delete(wave_frameset);
860  cpl_frameset_delete(gain_frameset);
861 
862 
863  gravi_data_delete(preproc_data);
864 
865  cpl_free(preproc_data);
866  cpl_free(calib_datas);
867  gravi_data_delete(profile_map);
868  gravi_data_delete(dark_map);
869  gravi_data_delete(wave_map);
870  gravi_data_delete(badpix_map);
871  gravi_data_delete(p2vm_map);
872  cpl_frameset_delete(recipe_frameset);
873 
874  return (int)cpl_error_set_message(cpl_func, CPL_ERROR_ILLEGAL_OUTPUT,
875  "Error while computing the preproc data_map");
876  }
877 
878  /* Option save the proproc file */
879  if (cpl_parameter_get_bool(p)){
880  preproc_name = cpl_sprintf("gravi_preproc_%s", filename);
881  cpl_frameset * preproc_frame = cpl_frameset_new ();
882  cpl_frameset_insert (preproc_frame, cpl_frame_duplicate (frame));
883 
884  applist = gravi_propertylist_get_qc (gravi_data_get_propertylist
885  (preproc_data, GRAVI_PRIMARY_HDR_NAME_EXT));
886 
887 
888  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, PREPROC);
889 
890 // gravi_data_save(preproc_data[i], frameset, preproc_name, parlist,
891 // recipe_frameset, "gravi_single", applist);
892  if (gravi_data_save(preproc_data, frameset, preproc_name, parlist,
893  preproc_frame, frame, "gravi_single", applist)
894  != CPL_ERROR_NONE){
895  cpl_frameset_delete(darkcalib_frameset);
896  cpl_frameset_delete(dark_frameset);
897  cpl_frameset_delete(wave_frameset);
898  cpl_frameset_delete(gain_frameset);
899  cpl_free(calib_datas);
900  gravi_data_delete(profile_map);
901  gravi_data_delete(dark_map);
902  gravi_data_delete(wave_map);
903  gravi_data_delete(p2vm_map);
904  cpl_frameset_delete(p2vm_frameset);
905  cpl_frameset_delete(recipe_frameset);
906  cpl_propertylist_delete(applist);
907 
908  gravi_data_delete(preproc_data);
909 
910  cpl_free(preproc_name);
911 
912  return (int) cpl_error_set_message(cpl_func,
913  CPL_ERROR_ILLEGAL_OUTPUT, "Could not save the preproc data"
914  " on the output file");
915  }
916 
917  cpl_free(preproc_name);
918  cpl_propertylist_delete (applist);
919  cpl_frameset_delete (preproc_frame);
920  }
921 
922  /* Compute the flux and visibilities for each telescope and
923  * per acquisition */
924  cpl_msg_info (NULL, "reduce the p2vm of the file %s",
925  cpl_frame_get_filename(frame));
926  p2vm_reduce = gravi_p2vm_reduce (preproc_data, p2vm_map, "gravi_single");
927  gravi_data_delete(preproc_data);
928  if (cpl_error_get_code()){
929  gravi_data_delete(p2vm_reduce);
930  gravi_data_delete(preproc_data);
931  gravi_data_delete(p2vm_map);
932  cpl_frameset_delete(p2vm_frameset);
933  cpl_frameset_delete(recipe_frameset);
934 
935 
936  return (int) cpl_error_set_message(cpl_func,
937  CPL_ERROR_ILLEGAL_OUTPUT, "error while reducing the p2vm");
938  }
939  char * p2vm_name = cpl_sprintf("single_reduced_%s",
940  cpl_frame_get_filename(frame));
941 
942  gravi_data_save_data(p2vm_reduce, p2vm_name, CPL_IO_CREATE);
943  cpl_free(p2vm_name);
944 
945  /* Compute the complex visibilty, squared visibity and the cloture phases */
946  /* To compute the visibility from a set of observable */
947 
948  /* by inverting the V2PM to compute the P2VM. The matrix inversion is
949  * done using the singular value decomposition method */
950  cpl_msg_info (NULL, "reduce the visibilities");
951  gravi_vis_reduce (oi_vis, p2vm_reduce, nb_vis, parlist, "gravi_single");
952 
953  nb_vis++;
954 
955  gravi_data_delete(p2vm_reduce);
956  }
957 
958  cpl_free(calib_datas);
959  gravi_data_delete(profile_map);
960  gravi_data_delete(dark_map);
961  gravi_data_delete(wave_map);
962  gravi_data_delete(badpix_map);
963  gravi_data_delete(p2vm_map);
964  cpl_parameterlist_delete (preproc_parlist);
965  cpl_frameset_delete(darkcalib_frameset);
966  cpl_frameset_delete(wave_frameset);
967  cpl_frameset_delete(dark_frameset);
968  cpl_frameset_delete(gain_frameset);
969  cpl_frameset_delete(badpix_frameset);
970  }
971  /* Initialization of the output name */
972 
973 
974  else{
975  /*test*/
976  recipe_frameset = cpl_frameset_duplicate (preproc_frameset);
977  nb_vis = cpl_frameset_get_size (recipe_frameset);
978  nb_frame = cpl_frameset_get_size(recipe_frameset);
979 
980  /* Load the p2vm data */
981  frame = cpl_frameset_get_position(p2vm_frameset, 0);
982  filename = cpl_frame_get_filename(frame);
983  cpl_msg_info (NULL, "This file %s is a p2vm file already "
984  "computed", filename);
985 
986  p2vm_map = gravi_data_load (filename);
987 
988  /* Compute the complex visibilty, squared visibity and the cloture phases */
989  /* To compute the visibility from a set of observable */
990 
991  /* by inverting the V2PM to compute the P2VM. The matrix inversion is
992  * done using the singular value decomposition method */
993 
994  oi_vis = gravi_data_new(0);
995 
996  for (i = 0; i < nb_frame; i++){
997  frame = cpl_frameset_get_position(preproc_frameset, i);
998  filename = cpl_frame_get_filename(frame);
999  cpl_msg_info (NULL, "Calibration of this file %s",
1000  filename);
1001  preproc_data = gravi_data_load(filename);
1002 
1003  cpl_msg_info (NULL, "reduce the p2vm of the file %s",
1004  cpl_frame_get_filename(frame));
1005  /* Compute the flux and visibilities for each telescope and
1006  * per acquisition */
1007  p2vm_reduce = gravi_p2vm_reduce (preproc_data, p2vm_map, "gravi_single");
1008  gravi_data_delete(preproc_data);
1009 
1010  if (cpl_error_get_code()){
1011 
1012  gravi_data_delete(p2vm_reduce);
1013  gravi_data_delete(preproc_data);
1014 
1015  gravi_data_delete(p2vm_map);
1016  cpl_frameset_delete(p2vm_frameset);
1017  cpl_frameset_delete(recipe_frameset);
1018 
1019 
1020  return (int) cpl_error_set_message(cpl_func,
1021  CPL_ERROR_ILLEGAL_OUTPUT, "error while reducing the p2vm");
1022  }
1023  char * p2vm_name = cpl_sprintf("single_reduced_%s",
1024  cpl_frame_get_filename(frame));
1025 
1026  gravi_data_save_data(p2vm_reduce, p2vm_name, CPL_IO_CREATE);
1027  cpl_free(p2vm_name);
1028 
1029 
1030  cpl_msg_info (NULL, "reduce the visibilities");
1031  gravi_vis_reduce (oi_vis, p2vm_reduce, i, parlist, "gravi_single");
1032 
1033  gravi_data_delete(p2vm_reduce);
1034 
1035  }
1036  gravi_data_delete(p2vm_map);
1037  }
1038 
1039 // char tag[20] = "VIS";
1040 
1041 
1042 
1043 // gravi_data_save_data(oi_vis, vis_name, CPL_IO_CREATE);
1044  /* Save the out put data file */
1045  frame = cpl_frameset_get_position (recipe_frameset, 0);
1046  frame_tag = cpl_frame_get_tag (frame);
1047  primary_hdr = gravi_data_get_propertylist(oi_vis,
1048  GRAVI_PRIMARY_HDR_NAME_EXT);
1049  applist = gravi_propertylist_get_qc (primary_hdr);
1050 
1051  cpl_propertylist_append_int(applist, GRAVI_NIGHT_OBS,
1052  cpl_propertylist_get_int (primary_hdr, GRAVI_NIGHT_OBS));
1053  cpl_propertylist_append_string(applist, "DATE-OBS",
1054  cpl_propertylist_get_string (primary_hdr, "DATE-OBS"));
1055 
1056  if ((strcmp(frame_tag, GRAVI_SINGLE_CALIB) == 0)) {
1057  filename = cpl_frame_get_filename(
1058  cpl_frameset_get_position (preproc_frameset, 0));
1059  vis_name = "vis_single_calib.fits";
1060  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_CALIB);
1061  }
1062  else if ((strcmp(frame_tag, GRAVI_SINGLE_SCIENCE) == 0)) {
1063  filename = cpl_frame_get_filename(
1064  cpl_frameset_get_position (preproc_frameset, 0));
1065  vis_name = "vis_single_science.fits";
1066  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_SCIENCE);
1067  }
1068  else {
1069  vis_name = "vis_single_science.fits";
1070  cpl_propertylist_append_string(applist, CPL_DFS_PRO_CATG, VIS_SCIENCE);
1071  }
1072 
1073  if (gravi_data_save(oi_vis, frameset, vis_name, parlist,
1074  recipe_frameset, frame, "gravi_single", applist)
1075  != CPL_ERROR_NONE){
1076 // gravi_data_delete(p2vm_map);
1077  gravi_data_delete(oi_vis);
1078  cpl_frameset_delete(p2vm_frameset);
1079  cpl_frameset_delete(recipe_frameset);
1080  cpl_propertylist_delete(applist);
1081 
1082  return (int) cpl_error_set_message(cpl_func,
1083  CPL_ERROR_ILLEGAL_OUTPUT, "Could not save the oi_vis"
1084  " on the output file");
1085  }
1086 
1087  cpl_propertylist_delete(applist);
1088 // if (gravi_save_qc (oi_vis, frameset, recipe_frameset,
1089 // vis_name, parlist, tag, "gravi_single") != CPL_ERROR_NONE){
1090 // printf("wayylouzzzzzzzze\n");
1091 // cpl_frameset_delete(wave_frameset);
1092 // cpl_frameset_delete(dark_frameset);
1093 // cpl_frameset_delete(gain_frameset);
1094 // gravi_data_delete(profile_map);
1095 // for (i = 0; i < nb_vis; i++){
1096 // gravi_data_delete(preproc_data[i]);
1097 // }
1098 // cpl_free(preproc_data);
1099 // for (i = 0; i < j; j++){
1100 // gravi_data_delete(p2vm_reduce[i]);
1101 // }
1102 // cpl_free(p2vm_reduce);
1103 // cpl_free(calib_datas);
1104 // gravi_data_delete(dark_map);
1105 // gravi_data_delete(wave_map);
1106 // gravi_data_delete(p2vm_map);
1107 // gravi_data_delete(oi_vis);
1108 // cpl_frameset_delete(p2vm_frameset);
1109 // gravi_data_delete(badpix_map);
1110 // cpl_frameset_delete(recipe_frameset);
1111 //
1112 // return (int)cpl_error_set_message(cpl_func, CPL_ERROR_ILLEGAL_OUTPUT,
1113 // "Error while save of the visibility file");
1114 // }
1115  gravi_data_delete(oi_vis);
1116 
1117 
1118  /* Deallocation of all variables */
1119 
1120  cpl_frameset_delete(p2vm_frameset);
1121  cpl_frameset_delete(recipe_frameset);
1122  cpl_frameset_delete(preproc_frameset);
1123 
1124 
1125  return (int)cpl_error_get_code();
1126 }
1127 
1128