67 const char * data,
const char *flag);
69 const char * data,
const char *flag);
72cpl_size
gravi_get_row_in_cat (cpl_table * diam_table,
double ra,
double dec,
double *separation);
88 cpl_ensure (array, CPL_ERROR_NULL_INPUT, -1);
91 if ( cpl_array_get_min ( array ) > 0.0 )
return 0;
93 cpl_size indx, size = cpl_array_get_size (array);
97 for ( indx = 0 ; indx < size ; indx ++ ) {
98 if ( cpl_array_get (array, indx, &nv) < 0.0 ) {
99 cpl_array_set_invalid (array, indx);
118 const char * data,
const char *flag)
121 cpl_ensure_code (oi_table, CPL_ERROR_NULL_INPUT);
122 cpl_ensure_code (data, CPL_ERROR_NULL_INPUT);
123 cpl_ensure_code (flag, CPL_ERROR_ILLEGAL_OUTPUT);
127 cpl_size row, nrow = cpl_table_get_nrow (oi_table);
128 cpl_array ** pdata = cpl_table_get_data_array (oi_table, data);
129 cpl_array ** pflag = cpl_table_get_data_array (oi_table, flag);
130 cpl_size indx, size = cpl_array_get_size (pdata[0]);
135 for ( row = 0 ; row < nrow ; row ++ ) {
136 for ( indx = 0 ; indx < size ; indx ++ ) {
137 if ( cpl_array_get (pdata[row], indx, &nv) < 0.0) {
138 cpl_array_set (pflag[row], indx, cpl_array_get (pflag[row], indx, &nv) + 1 );
144 return CPL_ERROR_NONE;
158 const char * data,
const char *flag)
161 cpl_ensure_code (oi_table, CPL_ERROR_NULL_INPUT);
162 cpl_ensure_code (data, CPL_ERROR_NULL_INPUT);
163 cpl_ensure_code (flag, CPL_ERROR_ILLEGAL_OUTPUT);
167 cpl_size row, nrow = cpl_table_get_nrow (oi_table);
168 cpl_array ** pdata = cpl_table_get_data_array(oi_table, data);
169 cpl_array ** pflag = cpl_table_get_data_array(oi_table, flag);
170 cpl_size indx, size = cpl_array_get_size (pdata[0]);
175 for ( row = 0 ; row < nrow ; row ++ ) {
176 if (pdata[row]==NULL)
continue;
177 for ( indx = 0 ; indx < size ; indx ++ ) {
178 if ( !cpl_array_is_valid (pdata[row], indx) ) {
179 cpl_array_set (pflag[row], indx, cpl_array_get (pflag[row], indx, &nv) + 1 );
185 return CPL_ERROR_NONE;
204 cpl_ensure (data, CPL_ERROR_NULL_INPUT, NULL);
210 char* setupstring = cpl_sprintf(
"%10s %s %s %.2fms %.2fs",
217 cpl_msg_debug (cpl_func,
"Get setup string: %s", setupstring);
252 const char * extName,
253 const char * insName,
254 const char * ampName,
255 const char * ampErrName,
256 int nbase,
double delta_t)
259 cpl_ensure_code (science, CPL_ERROR_NULL_INPUT);
260 cpl_ensure_code (insName, CPL_ERROR_NULL_INPUT);
261 cpl_ensure_code (extName, CPL_ERROR_NULL_INPUT);
262 cpl_ensure_code (ampName, CPL_ERROR_NULL_INPUT);
263 cpl_ensure_code (ampErrName, CPL_ERROR_NULL_INPUT);
264 cpl_ensure_code (used_tf_data, CPL_ERROR_NULL_INPUT);
265 cpl_ensure_code (num_tf_data>0, CPL_ERROR_ILLEGAL_INPUT);
267 cpl_msg_debug (cpl_func,
"%s %s amp=%s ampErr=%s nbase=%i",extName,insName,ampName,ampErrName,nbase);
269 int i, row_sc, row_cal, nv=0;
276 int nrow_sc = cpl_table_get_nrow (sci_table);
277 int nwave = cpl_table_get_column_depth (sci_table, ampName);
282 for (row_sc = 0; row_sc < nrow_sc; row_sc ++){
283 double time_sci = cpl_table_get_double (sci_table,
"TIME", row_sc, &nv);
284 int base = row_sc % nbase;
287 cpl_array * tf_mean = cpl_array_new (nwave, CPL_TYPE_DOUBLE);
288 cpl_array_fill_window (tf_mean, 0, nwave, 0.0);
293 for (i = 0; i < num_tf_data; i++){
296 int nrow_tf = cpl_table_get_nrow (tf_table);
299 for (row_cal = base; row_cal < nrow_tf; row_cal += nbase){
300 double time_tf = cpl_table_get_double (tf_table,
"TIME", row_cal, &nv);
303 cpl_array * tf_data = cpl_array_duplicate ( cpl_table_get_array (tf_table, ampName, row_cal) );
304 cpl_array * tf_err = cpl_array_duplicate ( cpl_table_get_array (tf_table, ampErrName, row_cal) );
309 double sigma = cpl_array_get_median (tf_err);
310 sigma = CPL_MAX (sigma, 0.02*fabs(cpl_array_get_median (tf_data)));
311 sigma = CPL_MAX (sigma, 1e-10);
314 double coeff = exp (-2 * fabs (time_sci - time_tf) / delta_t) / (sigma*sigma);
315 cpl_array_multiply_scalar (tf_data, coeff);
316 cpl_array_add (tf_mean, tf_data);
319 cpl_array_delete (tf_data);
320 cpl_array_delete (tf_err);
328 cpl_array_divide_scalar (tf_mean, norm);
335 if (sci_tf_table) cpl_table_set_array (sci_tf_table, ampName, row_sc, tf_mean);
338 cpl_array_divide (cpl_table_get_data_array (sci_table, ampName)[row_sc], tf_mean);
340 cpl_array_divide (cpl_table_get_data_array (sci_table, ampErrName)[row_sc], tf_mean);
343 cpl_array_delete (tf_mean);
359 return CPL_ERROR_NONE;
389 const char * extName,
390 const char * insName,
391 const char * phiName,
392 const char * phiErrName,
393 int nbase,
double delta_t)
396 cpl_ensure_code (science, CPL_ERROR_NULL_INPUT);
397 cpl_ensure_code (insName, CPL_ERROR_NULL_INPUT);
398 cpl_ensure_code (extName, CPL_ERROR_NULL_INPUT);
399 cpl_ensure_code (phiName, CPL_ERROR_NULL_INPUT);
400 cpl_ensure_code (phiErrName, CPL_ERROR_NULL_INPUT);
401 cpl_ensure_code (used_tf_data, CPL_ERROR_NULL_INPUT);
402 cpl_ensure_code (num_tf_data>0, CPL_ERROR_ILLEGAL_INPUT);
404 cpl_msg_debug(cpl_func,
"%s %s phi=%s phiErr=%s nbase=%i",extName,insName,phiName,phiErrName,nbase);
406 int i, row_sc, row_cal, nv=0;
413 int nrow_sc = cpl_table_get_nrow (sci_table);
414 int nwave = cpl_table_get_column_depth (sci_table, phiName);
419 for (row_sc = 0; row_sc < nrow_sc; row_sc ++){
420 double time_sci = cpl_table_get_double (sci_table,
"TIME", row_sc, &nv);
421 int base = row_sc % nbase;
424 cpl_array * tf_mean = cpl_array_new (nwave, CPL_TYPE_DOUBLE_COMPLEX);
425 cpl_array_fill_window_complex (tf_mean, 0, nwave, 0.0 + 0.0*I);
429 for (i = 0; i < num_tf_data; i++){
432 int nrow_tf = cpl_table_get_nrow (tf_table);
435 for (row_cal = base; row_cal < nrow_tf; row_cal += nbase){
436 double time_tf = cpl_table_get_double (tf_table,
"TIME", row_cal, &nv);
439 cpl_array * tf_data =
gravi_array_cexp ( I*CPL_MATH_PI/180.0, cpl_table_get_array (tf_table, phiName, row_cal) );
440 cpl_array * tf_err = cpl_array_duplicate (cpl_table_get_array (tf_table, phiErrName, row_cal) );
445 double sigma = CPL_MAX (cpl_array_get_median (tf_err), 0.2);
448 double coeff = exp (-2 * fabs (time_sci - time_tf) / delta_t) / (sigma*sigma);
449 cpl_array_multiply_scalar (tf_data, coeff);
450 cpl_array_add (tf_mean, tf_data);
452 cpl_array_delete (tf_data);
453 cpl_array_delete (tf_err);
462 cpl_array * sci_data =
gravi_array_cexp ( I*CPL_MATH_PI/180.0, cpl_table_get_array (sci_table, phiName, row_sc) );
463 cpl_array_divide (sci_data, tf_mean);
464 cpl_array_arg (sci_data);
465 cpl_array_multiply_scalar (sci_data, 180./CPL_MATH_PI);
466 cpl_table_set_array (sci_table, phiName, row_sc, sci_data);
471 cpl_array_arg (tf_mean);
472 cpl_array_multiply_scalar (tf_mean, 180./CPL_MATH_PI);
473 if (sci_tf_table) cpl_table_set_array (sci_tf_table, phiName, row_sc, tf_mean );
478 cpl_array_delete (tf_mean);
479 cpl_array_delete (sci_data);
494 return CPL_ERROR_NONE;
523 const cpl_parameterlist * parlist)
526 cpl_ensure (vis_data, CPL_ERROR_NULL_INPUT, NULL);
527 cpl_ensure (tf_data, CPL_ERROR_NULL_INPUT, NULL);
528 cpl_ensure (parlist, CPL_ERROR_NULL_INPUT, NULL);
529 cpl_ensure (num_tf>0, CPL_ERROR_NULL_INPUT, NULL);
533 char * setup_science, * setup_tf;
536 cpl_ensure( (vis_data != NULL) && (tf_data != NULL), CPL_ERROR_NULL_INPUT, NULL );
539 cpl_ensure(phi_tf_data != NULL, CPL_ERROR_NULL_INPUT, NULL);
548 cpl_msg_info (cpl_func,
"Number of possible TF: %i", num_tf);
552 cpl_msg_info (cpl_func,
"Delta time to interpolate the TF: %f s (%f h)", delta_t, delta_t/3600.0);
555 cpl_msg_info (cpl_func,
"Force calibration of the SCI by CALs: %s", force_calib?
"T":
"F");
559 cpl_msg_info (cpl_func,
"Setup of SCIENCE: %s", setup_science);
561 CPLCHECK_NUL(
"Cannot build the setup string for SCIENCE");
565 for (i = 0; i < num_tf; i++){
571 if (!(strcmp (setup_tf, setup_science )) ) {
573 used_tf_data[num_used_tf] = tf_data[i];
575 cpl_msg_info (cpl_func,
"Setup of TF file: %s -> keep", setup_tf);
576 }
else if (force_calib) {
578 used_tf_data[num_used_tf] = tf_data[i];
580 cpl_msg_info (cpl_func,
"Setup of TF file: %s -> keep (force_calib)", setup_tf);
583 cpl_msg_info (cpl_func,
"Setup of TF file: %s -> discard", setup_tf);
589 if (num_used_tf == 0) {
590 cpl_free (used_tf_data);
591 cpl_error_set_message (cpl_func, CPL_ERROR_NULL_INPUT,
"No calib file with the same keywords");
600 if (!force_calib && (strcmp (setup_tf, setup_science )) ) {
601 cpl_error_set_message (cpl_func, CPL_ERROR_NULL_INPUT,
"Visphi calib file does not have the same keywords");
615 int type_data, ntype_data = 2;
616 for (type_data = 0; type_data < ntype_data ; type_data ++) {
620 for ( pol= 0 ; pol < npol ; pol++ ) {
626 "VIS2DATA",
"VIS2ERR", 6, delta_t);
634 "VISAMP",
"VISAMPERR", 6, delta_t);
642 "T3AMP",
"T3AMPERR", 4, delta_t);
651 "VISPHI",
"VISPHIERR", 6, delta_t);
656 "VISPHI",
"VISPHIERR", 6, delta_t);
664 "T3PHI",
"T3PHIERR", 4, delta_t);
671 if(cpl_table_has_column(sci_table,
"FLUXDATA"))
672 fluxcolumn = cpl_sprintf(
"%s",
"FLUXDATA");
674 fluxcolumn = cpl_sprintf(
"%s",
"FLUX");
676 cpl_msg_info(cpl_func,
"Calibrating flux column %s", fluxcolumn);
680 fluxcolumn,
"FLUXERR", 4, delta_t);
681 cpl_free(fluxcolumn);
691 cpl_free (used_tf_data);
692 cpl_free (setup_science);
713 if (lbd <=0)
return 0.0;
714 double x = CPL_MATH_PI * uv / lbd * (diam * 1e-3 / 3600 / 180 * CPL_MATH_PI);
715 return ( (x<=0) ? 1.0 : 2.0 * j1 (x) / x );
734 cpl_ensure (diam_table, CPL_ERROR_NULL_INPUT, -1);
737 const char ** p_ra = cpl_table_get_data_string_const (diam_table,
"RAJ2000");
738 const char ** p_dec = cpl_table_get_data_string_const (diam_table,
"DEJ2000");
739 cpl_size nrow = cpl_table_get_nrow (diam_table);
741 cpl_ensure (p_ra, CPL_ERROR_ILLEGAL_INPUT, -1);
742 cpl_ensure (p_dec, CPL_ERROR_ILLEGAL_INPUT, -1);
743 cpl_ensure (nrow>0, CPL_ERROR_ILLEGAL_INPUT, -1);
750 for (cpl_size row = 0; row < nrow ; row ++) {
758 double dist = acos ( sin (c_dec) * sin (dec) + cos (c_dec) * cos (dec) * cos (ra - c_ra) );
768 if (separation != NULL) *separation = dis0 / CPL_MATH_PI * 180 * 3600;
789 cpl_ensure (vis_data, CPL_ERROR_NULL_INPUT, NULL);
792 cpl_propertylist *
plist;
793 cpl_table * tf_vistable;
794 cpl_table * tf_vis2table;
795 cpl_table * oi_wavelength;
797 double diameter, diameter_err, lambda, r_vis, r_vis2, tf_v, tf_v2;
805 for (type_data = 0; type_data < 2; type_data ++) {
806 cpl_table * diam_table;
808 cpl_size row_cat = -1;
809 char * diameter_source = NULL;
822 if ( (row_cat > -1) && (sep < 5.0) ) {
823 cpl_msg_info (cpl_func,
"Find match in DIAMETER_CAT (best sep=%.2f arcsec)", sep);
824 diameter = cpl_table_get (diam_table,
"UDDK", row_cat, NULL);
825 diameter_err = cpl_table_get (diam_table,
"e_LDD", row_cat, NULL);
826 diameter_source = cpl_sprintf(
"%s",
"CAT");
828 cpl_msg_warning (cpl_func,
"No match in DIAMETER_CAT (best sep=%.2f arcsec), use the HEADER value instead", sep);
831 diameter_source = cpl_sprintf(
"%s",
"HDR");
835 cpl_msg_info (cpl_func,
"No DIAMETER_CAT, use the HEADER value if any");
838 diameter_source = cpl_sprintf(
"%s",
"HDR");
843 if ( diameter>=4.0 ) {
844 cpl_msg_warning (cpl_func,
"Diameter for %s target: %.3f (+-%.3f) mas. Expected ?",
GRAVI_TYPE(type_data), diameter, diameter_err);
847 cpl_msg_warning (cpl_func,
"Diameter for %s target: 0.0 mas. Probably wrong, check parameter value",
GRAVI_TYPE(type_data));
850 cpl_msg_info(cpl_func,
"Diameter for %s target: %.3f (+-%.3f) mas",
GRAVI_TYPE(type_data), diameter, diameter_err);
855 for (pol = 0; pol < npol; pol ++)
862 if ((tf_vistable == NULL) || (tf_vis2table == NULL) || (oi_wavelength == NULL)){
864 cpl_error_set_message (cpl_func, CPL_ERROR_NULL_INPUT,
865 "Missing OI_VIS or OI_VIS2 or OI_WAVELENGTH");
873 int nrow = cpl_table_get_nrow (tf_vistable);
874 int nwave = cpl_table_get_column_depth (tf_vis2table,
"VIS2DATA" );
875 for (cpl_size row = 0; row < nrow; row ++){
878 r_vis = sqrt (pow (cpl_table_get_double (tf_vistable,
"UCOORD", row, &nv), 2) +
879 pow (cpl_table_get_double (tf_vistable,
"VCOORD", row, &nv), 2));
881 r_vis2 = sqrt (pow (cpl_table_get_double (tf_vis2table,
"UCOORD", row, &nv), 2) +
882 pow (cpl_table_get_double (tf_vis2table,
"VCOORD", row, &nv), 2));
887 cpl_array * model_vis = cpl_array_new (nwave, CPL_TYPE_DOUBLE);
888 cpl_array * model_vis2 = cpl_array_new (nwave, CPL_TYPE_DOUBLE);
891 for (cpl_size wave = 0; wave < nwave ; wave ++){
893 lambda = cpl_table_get (oi_wavelength,
"EFF_WAVE", wave, &nv);
898 cpl_array_set_double (model_vis, wave, tf_v);
899 cpl_array_set_double (model_vis2, wave, tf_v2);
904 cpl_array * tf_vis = cpl_table_get_data_array (tf_vistable,
"VISAMP")[row];
905 cpl_array * tf_visErr = cpl_table_get_data_array (tf_vistable,
"VISAMPERR")[row];
906 cpl_array_divide (tf_vis, model_vis);
907 cpl_array_divide (tf_visErr, model_vis);
910 cpl_array * tf_vis2 = cpl_table_get_data_array (tf_vis2table,
"VIS2DATA")[row];
911 cpl_array * tf_vis2Err = cpl_table_get_data_array (tf_vis2table,
"VIS2ERR")[row];
912 cpl_array_divide (tf_vis2, model_vis2);
913 cpl_array_divide (tf_vis2Err, model_vis2);
936 cpl_propertylist_set_comment (
plist, qc_name,
"TF rel.err from diameter at 2.2um");
941 cpl_propertylist_set_comment (
plist, qc_name,
"TF. VIS median over lbd.");
946 cpl_propertylist_set_comment (
plist, qc_name,
"TF. VIS2 median over lbd.");
950 cpl_array_delete (model_vis2);
951 cpl_array_delete (model_vis);
958 sprintf (qc_name,
"ESO QC TF CAL_%s DIAM",
GRAVI_TYPE(type_data));
960 sprintf (qc_name,
"ESO QC TF CAL_%s DIAM ERR",
GRAVI_TYPE(type_data));
962 sprintf (qc_name,
"ESO QC TF CAL_%s DIAM SOURCE",
GRAVI_TYPE(type_data));
963 cpl_propertylist_update_string (
plist, qc_name, diameter_source);
964 cpl_free(diameter_source);
968 if ( (row_cat > -1) && (sep < 5.0) ) {
972 const char sta = *cpl_table_get_string (oi_array,
"TEL_NAME", 0);
973 if (*cpl_table_get_string (oi_array,
"TEL_NAME", 1) != sta ||
974 *cpl_table_get_string (oi_array,
"TEL_NAME", 2) != sta ||
975 *cpl_table_get_string (oi_array,
"TEL_NAME", 3) != sta)
977 cpl_error_set_message (cpl_func, CPL_ERROR_ILLEGAL_INPUT,
"AT/UT mode is not supported");
983 if (sta==
'U') { diam = 8.0;
cpl_msg_info (cpl_func,
"Use UTs"); }
984 else if (sta==
'A') { diam = 1.8;
cpl_msg_info (cpl_func,
"Use ATs"); }
985 else { diam = 1.0; cpl_msg_warning (cpl_func,
"Cannot find the diameter of telescope"); }
988 double Kmag = cpl_table_get (diam_table,
"Kmag", row_cat, NULL);
989 double flux0 = 1.71173e+09 * pow (10.0, -Kmag/2.5) * (CPL_MATH_PI * pow (diam/2, 2));
990 cpl_msg_info (cpl_func,
"Use Kmag=%.2f for QC.TRANS", Kmag);
993 for (
int tel = 0; tel < 4; tel++) {
997 sprintf (qc_name,
"ESO QC FLUXRATE_%s%i_P1 SUM",
GRAVI_TYPE(type_data),tel+1);
999 sprintf (qc_name,
"ESO QC FLUXRATE_%s%i_P2 SUM",
GRAVI_TYPE(type_data),tel+1);
1004 sprintf (qc_name,
"ESO QC TF TRANS_%s%i",
GRAVI_TYPE(type_data),tel+1);
1006 cpl_propertylist_set_comment (
plist, qc_name,
"[%] Total transmission");
1007 cpl_msg_info (cpl_func,
"%s = %.2f%%", qc_name, flux / flux0 * 100.0);
1010 sprintf (qc_name,
"ESO QC TF CAL_%s DIAM CATNAME",
GRAVI_TYPE(type_data));
1011 cpl_propertylist_update_string (
plist, qc_name, cpl_table_get_string(diam_table,
"Name", row_cat));
1038 cpl_ensure_code (oi_vis, CPL_ERROR_NULL_INPUT);
1040 cpl_errorstate prestate = cpl_errorstate_get();
1042 cpl_msg_info (cpl_func,
"Compute the QC TF parameters");
1048 if ( cpl_error_get_code() || oi_tf == NULL) {
1049 cpl_msg_warning (cpl_func,
"Cannot compute the QC TF parameters for this observation... continue.");
1050 cpl_errorstate_set (prestate);
1051 return CPL_ERROR_NONE;
1061 return CPL_ERROR_NONE;
1073 cpl_ensure (vis_calib, CPL_ERROR_NULL_INPUT, NULL);
1074 cpl_ensure (num_calib>0, CPL_ERROR_ILLEGAL_INPUT, NULL);
typedefCPL_BEGIN_DECLS struct _gravi_data_ gravi_data
#define gravi_data_is_internal(data)
#define gravi_data_get_header(data)
#define gravi_data_get_oi_vis2(data, type, pol, npol)
#define gravi_data_get_oi_wave(data, type, pol, npol)
#define gravi_data_get_oi_vis(data, type, pol, npol)
cpl_msg_debug(cpl_func, "Spectra has <50 pixels -> don't flat")
cpl_msg_info(cpl_func, "Compute WAVE_SCAN for %s", GRAVI_TYPE(type_data))
#define GRAVI_OI_VIS2_EXT
#define GRAVI_INSNAME(type, pol, npol)
#define GRAVI_OI_ARRAY_EXT
#define GRAVI_OI_FLUX_EXT
#define gravi_pfits_get_type_raep(plist, type)
#define gravi_pfits_get_type_decep(plist, type)
#define gravi_msg_function_exit(flag)
#define FREE(function, variable)
#define CPLCHECK_NUL(msg)
#define gravi_msg_function_start(flag)
#define CPLCHECK_MSG(msg)
cpl_array * gravi_array_cexp(double complex factor, const cpl_array *input)
Compute the complex exponention of an array: cexp (factor * input)
gravi_data * gravi_data_duplicate(const gravi_data *self)
Create a copy of the gravi data.
cpl_table * gravi_data_get_table(gravi_data *self, const char *extname)
Return a pointer on a table extension by its EXTNAME.
cpl_table * gravi_data_get_oi_table(gravi_data *self, const char *extname, const char *insname)
Get an OI_FITS table from EXTNAME and INSNAME.
cpl_table * gravi_data_get_table_x(gravi_data *self, int i)
Get the table of an extension by position.
void gravi_data_delete(gravi_data *self)
Delete a gravi data.
int gravi_param_get_bool(const cpl_parameterlist *parlist, const char *name)
double gravi_param_get_double(const cpl_parameterlist *parlist, const char *name)
int gravi_pfits_get_pola_num(const cpl_propertylist *plist, int type_data)
cpl_error_code gravi_pfits_update_double(cpl_propertylist *plist, const char *full_name, double value)
const char * gravi_pfits_get_pola_mode(const cpl_propertylist *plist, int type_data)
const char * gravi_pfits_get_spec_res(const cpl_propertylist *plist)
double gravi_pfits_get_dit_sc(const cpl_propertylist *plist)
double gravi_dec_to_rad(const char *stri)
double gravi_pfits_get_dit_ft(const cpl_propertylist *plist)
double gravi_pfits_get_double_silentdefault(const cpl_propertylist *plist, const char *name, double def)
double gravi_pfits_get_diameter(const cpl_propertylist *plist, int type_data)
double gravi_ra_to_rad(const char *stri)
char * gravi_calib_setupstring(gravi_data *data)
Build a unique setup string per calibratable setup.
gravi_data * gravi_calibrate_vis(gravi_data *vis_data, gravi_data **tf_data, int num_tf, gravi_data *phi_tf_data, gravi_data *tf_science, const cpl_parameterlist *parlist)
Computes the calibrated visibility from science a single data and several previously evaluated instru...
cpl_error_code gravi_vis_flag_negative(cpl_table *oi_table, const char *data, const char *flag)
Flag negative element of an OIFITS table.
double gravi_visibility_UD(double uv, double diam, double lbd)
Compute the expected visibility from a UD model.
int gravi_array_set_invalid_negative(cpl_array *array)
Set invalid to all negative elements of an array.
gravi_data * gravi_compute_zp(gravi_data **vis_calib, int num_calib)
Compute the ZP data.
cpl_size gravi_get_row_in_cat(cpl_table *diam_table, double ra, double dec, double *separation)
Get the row in a cpl_table matching closest RAJ2000 and DEJ2000.
cpl_error_code gravi_apply_tf_amp(gravi_data *science, gravi_data *science_tf, gravi_data **used_tf_data, int num_tf_data, const char *extName, const char *insName, const char *ampName, const char *ampErrName, int nbase, double delta_t)
Interpolate the TF at the time of the science observation for an amplitude quantity.
cpl_error_code gravi_compute_tf_qc(gravi_data *oi_vis, gravi_data *diamcat_data)
Fill QC parameters related to transfer function.
cpl_error_code gravi_apply_tf_phi(gravi_data *science, gravi_data *science_tf, gravi_data **used_tf_data, int num_tf_data, const char *extName, const char *insName, const char *phiName, const char *phiErrName, int nbase, double delta_t)
Interpolate the TF at the time of the science observation for a phase quantity (in deg).
cpl_error_code gravi_vis_flag_invalid(cpl_table *oi_table, const char *data, const char *flag)
Flag invalid element of an OIFITS table.
gravi_data * gravi_compute_tf(gravi_data *vis_data, gravi_data *diamcat_data)
This function evaluates the transfer function from the observation of a reference star whose diameter...
char GRAVI_BASE_NAME[6][3]