/* * This file is part of the ESO SINFONI Pipeline * Copyright (C) 2004,2005 European Southern Observatory * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by * the Free Software Foundation; either version 2 of the License, or * (at your option) any later version. * * This program is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * GNU General Public License for more details. * * You should have received a copy of the GNU General Public License * along with this program; if not, write to the Free Software * Foundation, 51 Franklin St, Fifth Floor, Boston, MA 02111-1307 USA */ /*---------------------------------------------------------------------------- File name : sinfo_find_distortions.c Author : A. Modigliani Created on : Aug 12, 2004 Description : ---------------------------------------------------------------------------*/ #ifdef HAVE_CONFIG_H # include #endif #include /*---------------------------------------------------------------------------- Includes ---------------------------------------------------------------------------*/ #include "sinfo_new_find_distortions.h" #include "sinfo_pro_save.h" #include "sinfo_pro_types.h" #include "sinfo_finddist_ini_by_cpl.h" #include "sinfo_wave_calibration.h" #include "sinfo_cube_construct.h" #include "sinfo_absolute.h" #include "sinfo_distortion.h" #include "sinfo_utilities.h" #include "sinfo_recipes.h" #include "sinfo_utils_wrappers.h" #include "sinfo_error.h" #include "sinfo_globals.h" /*---------------------------------------------------------------------------- Defines ---------------------------------------------------------------------------*/ #define SINFO_ARC_SATURATION 100000 #define SINFO_ARC_MAX_WIDTH 64 /**@{*/ /** * @addtogroup sinfo_rec_mflat Master flat generation * * TBD */ /*---------------------------------------------------------------------------- Function Definitions ---------------------------------------------------------------------------*/ /*---------------------------------------------------------------------------- Function : sinfo_find_distortions() In : Out : Job : ---------------------------------------------------------------------------*/ static double new_compute_shift(double x, double y, double c00, double c10, double c01, double c11, double c20, double c02, double c21, double c12, double c30, double c03); int sinfo_new_find_distortions(const char* plugin_id, cpl_parameterlist* config, cpl_frameset* sof, cpl_frameset* set_fibre_ns) { finddist_config * cfg=NULL ; cpl_image * imonind=NULL ; cpl_image * impoly=NULL ; cpl_image * im=NULL ; cpl_image * map=NULL ; cpl_image * mask=NULL ; int* degx=NULL; int* degy=NULL; int* n_found_lines=NULL; int* sum_pointer=NULL; int** row_clean=NULL; int pdensity=0; int x_l=SIZEX/4*1; int x_u=SIZEX/4*3; int y_l=SIZEY/4*1; int y_u=SIZEY/4*3; int x_c=SIZEX/4*2; int y_c=SIZEY/4*2; const int pdx=4; const int pdy=4; int check = 0; int i = 0; int lx=0; //int ly=0; int sum=0; int n_lines=0; int fit=0; int n=0; int j=0; cpl_size coef_pow[2]; int arcs_window_size=0; float value=0; float newval=0; float total_dist=0; float shift=0; float* distances=NULL; float** wavelength_clean=NULL; float** slit_pos=NULL; float** acoefs=NULL; float* wave=NULL; float* intens=NULL; float* first =NULL; double* coef=NULL; double xshift=0.; double arcs_thres_factor=0; double arcs_min_arclen_factor=0; double pcf[pdx][pdy]; char tbl_name[FILE_NAME_SZ]; char key_name[FILE_NAME_SZ]; cpl_table* poly_tbl=NULL; cpl_table* tbl_spos=NULL; cpl_table* tbl_line_list=NULL; FitParams** par=NULL; cpl_frameset* stk=NULL; cpl_table* qclog_tbl=NULL; cpl_polynomial* distor_poly=NULL; cpl_apertures * arcs=NULL ; cpl_parameter* p=NULL; int smooth_rad=0; check_nomsg(p=cpl_parameterlist_find(config,"sinfoni.product.density")); check_nomsg(pdensity=cpl_parameter_get_int(p)); /* parse the file names and parameters to the finddist_config data structure cfg */ check_nomsg(stk=cpl_frameset_new()); cknull(cfg = sinfo_parse_cpl_input_finddist(config,sof,&stk), "could not parse CPL input!"); if(sinfo_is_fits_file (cfg->inFrame) != 1) { sinfo_msg_error("Input file %s is not FITS",cfg->inFrame); goto cleanup; } /* load the emission line frame--- */ check(im = cpl_image_load(cfg->inFrame,CPL_TYPE_FLOAT,0,0), "could not load arc image"); /* load the fake on - fake off frame--- */ check(imonind = cpl_image_load(cfg->nsFrame,CPL_TYPE_FLOAT,0,0), "could not load on-off fake image"); check(impoly = cpl_image_load(cfg->nsFrame,CPL_TYPE_FLOAT,0,0), "could not load on-off fake image"); /* load the fake on - fake off frame--- */ check(mask = cpl_image_load(cfg->mask,CPL_TYPE_FLOAT,0,0), "could not load mask image"); check(cpl_image_multiply (im,mask), "Failing to correct arc ima by bp mask"); check(cpl_image_multiply (imonind,mask), "Failing to correct on-off fake ima by bp mask"); check(cpl_image_multiply(impoly,mask), "Failing to correct on-off fake ima by bp mask"); sinfo_free_image(&mask); check_nomsg(lx = cpl_image_get_size_x(im)); //check_nomsg(ly = cpl_image_get_size_y(im)); /* TO BE CHENGED THE FOLLOWING INPUT */ /* open the line list and read the number of lines */ check(tbl_line_list = cpl_table_load(cfg->lineList,1,0), "problems loading table %s",cfg->lineList); check_nomsg(n_lines = cpl_table_get_nrow(tbl_line_list)); check_nomsg(wave=cpl_table_get_data_float(tbl_line_list,"wave")); cpl_type type=cpl_table_get_column_type(tbl_line_list,"int"); if(type==CPL_TYPE_INT) { check_nomsg(cpl_table_cast_column(tbl_line_list,"int", "fint",CPL_TYPE_FLOAT)); check_nomsg(intens = cpl_table_get_data_float(tbl_line_list,"fint")); } else { check_nomsg(intens = cpl_table_get_data_float(tbl_line_list,"int")); } /* #--------------------------------------------------------------------- #---------------------------FINDLINES--------------------------------- #--------------------------------------------------------------------- */ sinfo_msg("Find Lines"); /* do the wavelength calibration, that means: find the dispersion relation and parameterize its coefficients */ acoefs = sinfo_new_2Dfloatarray(cfg->nrDispCoefficients, lx); /* allocate memory */ n_found_lines = sinfo_new_intarray(lx); row_clean = sinfo_new_2Dintarray(lx, n_lines); wavelength_clean = sinfo_new_2Dfloatarray(lx, n_lines); sum_pointer = sinfo_new_intarray(1); /* find the emission lines in each image column */ sinfo_new_intarray_set_value(sum_pointer, 0, 0); ck0((check = sinfo_new_find_lines(im, wave, intens, n_lines, row_clean, wavelength_clean, cfg->guessBeginWavelength, cfg->guessDispersion1, cfg->guessDispersion2, cfg->mindiff, cfg->halfWidth, n_found_lines, cfg->sigma, sum_pointer )), "FindLines failed!"); sinfo_free_table(&tbl_line_list); /* #------------------------------------------------------------------------- #---------------------------WAVECALIB------------------------------------- #------------------------------------------------------------------------- */ sinfo_msg("Do wave calib"); sum = sinfo_new_intarray_get_value(sum_pointer,0); /* #allocate memory for the fit parameters */ cknull(par = sinfo_new_fit_params( sum ), "sinfo_new_fit_params failed!"); /* fit each line, make a polynomial fit and fit the resulting fit coefficients across the columns of the slitlet */ cknull(map = sinfo_new_wave_cal(im, par, acoefs, cfg->nslitlets, row_clean, wavelength_clean, n_found_lines, cfg->guessDispersion1, cfg->halfWidth, cfg->minAmplitude, cfg->maxResidual, cfg->fwhm, cfg->nrDispCoefficients, cfg->nrCoefCoefficients, cfg->sigmaFactor, cfg->pixeldist, cfg->pixel_tolerance), "sinfo_waveCal failed!"); /* here mem leak */ sinfo_msg("Check line positions"); shift = sinfo_new_check_line_positions (im, acoefs, cfg->nrDispCoefficients, cfg->guessDispersion1,par); if (FLAG == shift) { sinfo_msg_error ( "checkForLinePositions failed!"); goto cleanup; } sinfo_free_image(&map); /* free memory */ sinfo_new_destroy_2Dfloatarray ( &wavelength_clean, lx ); sinfo_new_destroy_2Dintarray (&row_clean, lx); sinfo_new_destroy_intarray(&n_found_lines ); sinfo_new_destroy_intarray(&sum_pointer ); sinfo_new_destroy_2Dfloatarray ( &acoefs, cfg->nrDispCoefficients ); /* sinfo_new_destroy_array (& wave ); sinfo_new_destroy_array (& intens ); */ /* #---------------------------------------------------------------------------- #-------------------------SLITFITS------------------------------------------- #---------------------------------------------------------------------------- #--fit the slitlet sinfo_edge positions if desired-- */ sinfo_msg("Find slit pos"); /*allocate memory for the slitlet position array */ slit_pos = sinfo_new_2Dfloatarray(32,2); fit = sinfo_new_fit_slits_boltz( im, par, slit_pos, cfg->boxLength, cfg->yBox, cfg->diffTol ); if (fit < 0) { sinfo_msg_error( "sinfo_fitSlitsBoltz failed" ); goto cleanup; } sinfo_new_destroy_fit_params(&par); check_nomsg(tbl_spos=cpl_table_new(32)); check_nomsg(cpl_table_new_column(tbl_spos,"pos1", CPL_TYPE_DOUBLE)); check_nomsg(cpl_table_new_column(tbl_spos,"pos2", CPL_TYPE_DOUBLE)); check_nomsg(cpl_table_set_column_format(tbl_spos,"pos1", "15.9f")); check_nomsg(cpl_table_set_column_format(tbl_spos,"pos2", "15.9f")); for (i =0; i< 32; i++) { check_nomsg(cpl_table_set_double(tbl_spos,"pos1",i, sinfo_new_array2D_get_value(slit_pos,i,0))); check_nomsg(cpl_table_set_double(tbl_spos,"pos2",i, sinfo_new_array2D_get_value(slit_pos,i,1))); } if(pdensity > 2) { strcpy(tbl_name,"out_slitlets_pos_predist.fits"); ck0(sinfo_pro_save_tbl(tbl_spos,set_fibre_ns,sof,tbl_name, PRO_SLITLETS_POS_PREDIST,NULL,plugin_id,config), "cannot save tbl %s", tbl_name); } sinfo_free_table(&tbl_spos); /* #--------------------------------------------------------- # do the north - south - test #--------------------------------------------------------- */ sinfo_msg("Do north south test"); /* sinfo_msg("cfg->nslits =%d\n", cfg->nslits); sinfo_msg("cfg->nshalfWidth =%d\n", cfg->nshalfWidth); sinfo_msg("cfg->nsfwhm=%f\n",cfg->nsfwhm ); sinfo_msg("cfg->minDiff=%f\n", cfg->minDiff); sinfo_msg("cfg->estimated_dist=%f\n", cfg->estimated_dist); sinfo_msg("cfg->devtol=%f\n", cfg->devtol); sinfo_msg("cfg->loPos=%d\n", cfg->loPos); sinfo_msg("cfg->hiPos=%d\n",cfg->hiPos); */ cknull(distances = sinfo_north_south_test(imonind, cfg->nslits, cfg->nshalfWidth, cfg->nsfwhm , cfg->minDiff, cfg->estimated_dist, cfg->devtol, cfg->loPos, cfg->hiPos), "north_south_test failed"); sinfo_free_image(&imonind); /* #--------------------------------------------------------- # get firstcol #--------------------------------------------------------- */ sinfo_msg("get first col"); first = sinfo_new_floatarray(61); total_dist=0; n=0; for (i=0; i< 31; i++) { total_dist=total_dist+sinfo_new_array_get_value(distances,i); for (j=0; j< 2; j++) { if (j == 0) { if (i != 30) { newval = sinfo_new_array2D_get_value(slit_pos,i+1,j) - total_dist; sinfo_new_array_set_value(first, newval, n); n = n+1; } } if (j == 1) { newval = sinfo_new_array2D_get_value(slit_pos,i,j) - total_dist; sinfo_new_array_set_value(first, newval, n); n = n+1; } } } value = sinfo_new_f_median(first,61); sinfo_msg("Firstcol =%f", value); sinfo_new_destroy_array(&first); sinfo_new_destroy_array(&distances); sinfo_new_destroy_2Dfloatarray ( &slit_pos, 32 ); /* #--------------------------------------------------------- # Determine distortions #--------------------------------------------------------- */ sinfo_msg("Determine distortions"); degx=cpl_calloc(8,sizeof(int)); degy=cpl_calloc(8,sizeof(int)); coef=cpl_calloc(8,sizeof(double)); check_nomsg(p=cpl_parameterlist_find(config, "sinfoni.distortion.arcs_thresh_factor")); check_nomsg(arcs_thres_factor=cpl_parameter_get_double(p)); check_nomsg(p=cpl_parameterlist_find(config, "sinfoni.distortion.arcs_min_arclen_factor")); check_nomsg(arcs_min_arclen_factor=cpl_parameter_get_double(p)); check_nomsg(p=cpl_parameterlist_find(config, "sinfoni.distortion.arcs_window_size")); check_nomsg(arcs_window_size=cpl_parameter_get_int(p)); check_nomsg(p=cpl_parameterlist_find(config, "sinfoni.distortion.smooth_rad")); check_nomsg(smooth_rad=cpl_parameter_get_int(p)); cknull(distor_poly=sinfo_distortion_estimate_new(impoly, 0, 0, cpl_image_get_size_x(impoly), cpl_image_get_size_y(impoly), FALSE, SINFO_ARC_SATURATION, SINFO_ARC_MAX_WIDTH, arcs_thres_factor, arcs_min_arclen_factor, arcs_window_size, smooth_rad, 3, (double)value, &arcs), "cannot estimage distortion") ; /*AMo: additional mem leaks */ sinfo_free_apertures(&arcs); coef_pow[0]=0; coef_pow[1]=0; check_nomsg(pcf[0][0]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",0,0,pcf[0][0]); /* pcf[0][0]+=value; */ sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",0,0,pcf[0][0]); check_nomsg(cpl_polynomial_set_coeff(distor_poly,coef_pow,pcf[0][0])); /* up to here ok */ /* To check mem leaks */ sinfo_free_image(&impoly); //sinfo_msg("Polynomial fit results: deg=%d", // cpl_polynomial_get_degree(distor_poly)); coef_pow[0]=0; coef_pow[1]=0; check_nomsg(pcf[0][0]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",0,0,pcf[0][0]); coef_pow[0]=1; coef_pow[1]=0; check_nomsg(pcf[1][0]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",1,0,pcf[1][0]); coef_pow[0]=0; coef_pow[1]=1; check_nomsg(pcf[0][1]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",0,1,pcf[0][1]); coef_pow[0]=1; coef_pow[1]=1; check_nomsg(pcf[1][1]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",1,1,pcf[1][1]); coef_pow[0]=2; coef_pow[1]=0; check_nomsg(pcf[2][0]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",2,0,pcf[2][0]); coef_pow[0]=0; coef_pow[1]=2; check_nomsg(pcf[0][2]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",0,2,pcf[0][2]); coef_pow[0]=2; coef_pow[1]=1; check_nomsg(pcf[2][1]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",2,1,pcf[2][1]); coef_pow[0]=1; coef_pow[1]=2; check_nomsg(pcf[1][2]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",1,2,pcf[1][2]); coef_pow[0]=3; coef_pow[1]=0; check_nomsg(pcf[3][0]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",3,0,pcf[3][0]); coef_pow[0]=0; coef_pow[1]=3; check_nomsg(pcf[0][3]=cpl_polynomial_get_coeff(distor_poly,coef_pow)); //sinfo_msg("Polynomial fit results: coeff[%d][%d]=%g",0,3,pcf[0][3]); sinfo_msg("Save results"); check_nomsg(poly_tbl=cpl_table_new(10)); check_nomsg(cpl_table_new_column(poly_tbl,"degx", CPL_TYPE_INT)); check_nomsg(cpl_table_new_column(poly_tbl,"degy", CPL_TYPE_INT)); check_nomsg(cpl_table_new_column(poly_tbl,"coeff", CPL_TYPE_DOUBLE)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",0,0)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",1,1)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",2,0)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",3,1)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",4,2)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",5,0)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",6,2)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",7,1)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",8,3)); check_nomsg(cpl_table_set_int(poly_tbl,"degx",9,0)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",0,0)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",1,0)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",2,1)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",3,1)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",4,0)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",5,2)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",6,1)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",7,2)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",8,0)); check_nomsg(cpl_table_set_int(poly_tbl,"degy",9,3)); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",0,pcf[0][0])); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",1,pcf[1][0])); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",2,pcf[0][1])); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",3,pcf[1][1])); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",4,pcf[2][0])); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",5,pcf[0][2])); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",6,pcf[2][1])); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",7,pcf[1][2])); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",8,pcf[3][0])); check_nomsg(cpl_table_set_double(poly_tbl,"coeff",9,pcf[0][3])); //sinfo_msg("Polynomial distortion coefficients \n"); //sinfo_msg("c= %g %g %g %g %g %g %g %g %g %g\n", // pcf[0][0],pcf[1][0],pcf[0][1],pcf[1][1], // pcf[2][0],pcf[0][2],pcf[2][1],pcf[1][2],pcf[3][0],pcf[0][3]); /* QC LOG */ check_nomsg(qclog_tbl = sinfo_qclog_init()); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",0,0); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[0][0], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",1,0); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[1][0], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",0,1); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[0][1], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",1,1); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[1][1], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",2,0); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[2][0], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",0,2); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[0][2], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",2,1); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[2][1], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",1,2); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[1][2], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",3,0); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[3][0], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s%d%d","QC COEFF",0,3); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,pcf[0][3], "Polynomial distortion coefficient","%g")); snprintf(key_name,MAX_NAME_SIZE-1,"%s","QC OFFSET"); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,value, "Polynomial distortion coefficient","%g")); xshift=new_compute_shift(x_c,y_c,pcf[0][0],pcf[1][0],pcf[0][1], pcf[1][1],pcf[2][0],pcf[0][2], pcf[2][1],pcf[1][2],pcf[3][0],pcf[0][3]); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC XSHIFT CC",xshift, "X shift in x_c,y_c","%g")); xshift=new_compute_shift(x_l,y_l,pcf[0][0],pcf[1][0],pcf[0][1], pcf[1][1],pcf[2][0],pcf[0][2], pcf[2][1],pcf[1][2],pcf[3][0],pcf[0][3]); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC XSHIFT LL",xshift, "X shift in x_l,y_l","%g")); xshift=new_compute_shift(x_l,y_u,pcf[0][0],pcf[1][0],pcf[0][1], pcf[1][1],pcf[2][0],pcf[0][2], pcf[2][1],pcf[1][2],pcf[3][0],pcf[0][3]); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC XSHIFT UL",xshift, "X shift in x_l,y_u","%g")); xshift=new_compute_shift(x_u,y_u,pcf[0][0],pcf[1][0],pcf[0][1], pcf[1][1],pcf[2][0],pcf[0][2], pcf[2][1],pcf[1][2],pcf[3][0],pcf[0][3]); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC XSHIFT UR",xshift, "X shift in x_u,y_u","%g")); xshift=new_compute_shift(x_u,y_l,pcf[0][0],pcf[1][0],pcf[0][1], pcf[1][1],pcf[2][0],pcf[0][2], pcf[2][1],pcf[1][2],pcf[3][0],pcf[0][3]); ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC XSHIFT LR",xshift, "X shift in x_u,y_l","%g")); ck0(sinfo_pro_save_tbl(poly_tbl,set_fibre_ns,sof,cfg->outName, PRO_DISTORTION,qclog_tbl,plugin_id,config), "cannot dump tbl %s", cfg->outName); sinfo_free_table(&poly_tbl); sinfo_free_table(&qclog_tbl); sinfo_free_polynomial(&distor_poly); sinfo_free_int(°x); sinfo_free_int(°y); sinfo_free_double(&coef); sinfo_free_image(&im); sinfo_free_frameset(&stk); sinfo_finddist_free (&cfg); return 0; cleanup: sinfo_free_table(&poly_tbl); sinfo_free_table(&qclog_tbl); sinfo_free_polynomial(&distor_poly); sinfo_free_int(°x); sinfo_free_int(°y); sinfo_free_double(&coef); sinfo_free_apertures(&arcs); /*if(wave != NULL) sinfo_new_destroy_array (& wave );*/ /*if(intens != NULL) sinfo_new_destroy_array (& intens );*/ if(first != NULL) sinfo_new_destroy_array(&first); sinfo_free_table(&tbl_spos); if(slit_pos != NULL) sinfo_new_destroy_2Dfloatarray (&slit_pos,32); if(distances != NULL) sinfo_new_destroy_array(&distances); if(par!=NULL) sinfo_new_destroy_fit_params(&par); if(n_found_lines != NULL) sinfo_new_destroy_intarray(&n_found_lines ); if(row_clean != NULL) sinfo_new_destroy_2Dintarray(&row_clean, lx); if(wavelength_clean != NULL) { sinfo_new_destroy_2Dfloatarray(&wavelength_clean,lx); } if(sum_pointer != NULL) sinfo_new_destroy_intarray(&sum_pointer); if(acoefs != NULL) { sinfo_new_destroy_2Dfloatarray(&acoefs, cfg->nrDispCoefficients ); } sinfo_free_table(&tbl_line_list); sinfo_free_image(&mask); sinfo_free_image(&impoly); sinfo_free_image(&imonind); sinfo_free_image(&map); sinfo_free_image(&im); sinfo_finddist_free (&cfg); sinfo_free_frameset(&stk); return -1; } static double new_compute_shift(double x, double y, double c00, double c10, double c01, double c11, double c20, double c02, double c21, double c12, double c30, double c03) { double x_shift=0; double shift=0; //sinfo_msg("c= %g %g %g %g %g %g %g %g %g %g\n", // c00,c10,c01,c11, // c20,c02,c21,c12,c30,c03); shift=c00+c10*x+c01*y+ c11*x*y+c20*x*x+c02*y*y+ c21*x*x*y+c12*x*y*y+c30*x*x*x+c03*y*y*y; x_shift=x-shift; return x_shift; } /**@}*/