find_distortions.c

00001 /*----------------------------------------------------------------------------
00002  
00003    File name    :       find_distortions.c
00004    Author       :   A. Modigliani
00005    Created on   :   Aug 12, 2004
00006    Description  : 
00007 
00008  ---------------------------------------------------------------------------*/
00009 
00010 /*----------------------------------------------------------------------------
00011                                 Includes
00012  ---------------------------------------------------------------------------*/
00013 #include "find_distortions.h"
00014 #include "sinfoni_pro_save.h"
00015 #include "finddist_ini_by_cpl.h"
00016 
00017 /*----------------------------------------------------------------------------
00018                                 Defines
00019  ---------------------------------------------------------------------------*/
00020 
00021 
00022 /*----------------------------------------------------------------------------
00023                              Function Definitions
00024  ---------------------------------------------------------------------------*/
00025 
00026 /*----------------------------------------------------------------------------
00027    Function     :       find_distortions()
00028    In           :       
00029    Out          :        
00030    Job          :
00031 
00032  ---------------------------------------------------------------------------*/
00033 int
00034 find_distortion_qclog(cpl_table* poly_tbl, cpl_table* qclog_tbl);
00035 double
00036 compute_shift(double x,
00037               double y,
00038               double c00,
00039               double c10,
00040               double c01,
00041               double c11,
00042               double c20,
00043               double c02,
00044               double c21,
00045               double c12);
00046 
00047 int find_distortions (const char* plugin_id, cpl_parameterlist* config, cpl_frameset* sof)
00048 {
00049     const char* _id = "find_distortions";
00050     finddist_config * cfg ;
00051     OneImage * imonind=NULL ;
00052     OneImage * impoly=NULL ;
00053     OneImage * im=NULL ;
00054     OneImage *  map=NULL ;
00055     OneImage *  mask=NULL ;
00056     int*            degx;
00057     int*            degy;
00058     double*         coef;
00059 
00060     char* poly_file=NULL;
00061     char        ch,temp[150];
00062     char        poly_string[150];
00063     int     count=0;
00064     FILE    *   poly_in;
00065 
00066      
00067 
00068     int x_l=SIZEX/4*1;
00069     int x_u=SIZEX/4*3;
00070     int y_l=SIZEY/4*1;
00071     int y_u=SIZEY/4*3;
00072     int x_c=SIZEX/4*2;
00073     int y_c=SIZEY/4*2;
00074 
00075   
00076     int check = 0;
00077     int i = 0;
00078     int lx=0;
00079     int ly=0;
00080     poly2d* poly=NULL;
00081 
00082     float value=0;
00083     float* distances=NULL;
00084     float** wavelength_clean;
00085     float** slit_pos;
00086     float** acoefs;
00087     float newval;
00088   int    dx0,dy0,dx1,dy1,dx2,dy2;
00089   float   c00,c10,c01,c11,c20,c02,c21,c12;
00090 
00091     float* wave=NULL;
00092     float* intens=NULL;
00093     float* first =NULL;
00094     cpl_table* poly_tbl=NULL;
00095     float total_dist=0;
00096     float shift=0;
00097     double xshift=0.;
00098 
00099     int sum=0;
00100     int n_lines=0;
00101     int fit=0;
00102     int n=0;
00103     int j=0;
00104 
00105     char* key_value=NULL;
00106     char* key_name=NULL;
00107     cpl_table* tbl_spos=NULL;
00108 
00109     cpl_table* tbl_line_list=NULL;
00110     FitParams** par;
00111 
00112     cpl_frameset* stk=NULL;
00113     cpl_table* qclog_tbl=NULL;
00114 
00115 
00116     int** row_clean;
00117     int* n_found_lines=NULL;
00118     int* sum_pointer=NULL;
00119     char tbl_name[FILE_NAME_SZ];
00120 
00121     /* parse the file names and parameters to the finddist_config 
00122        data structure cfg 
00123      */
00124     stk=cpl_frameset_new();
00125     cfg = parse_cpl_input_finddist(config,sof,&stk);
00126 
00127 
00128     if (cfg == NULL) {
00129       cpl_msg_error (_id,"could not parse CPL input!");
00130       cpl_frameset_delete(stk);
00131       return -1;
00132     }
00133 
00134     if(is_fits_file (cfg->inFrame) != 1) {
00135       cpl_msg_error(_id,"Input file %s is not FITS",cfg->inFrame);
00136       finddist_cfg_destroy (cfg);
00137       cpl_frameset_delete(stk);
00138       return -1;
00139     }
00140 
00141 
00142     /* load the emission line frame--- */
00143     im = load_image(cfg->inFrame);                
00144     if (im == NULL) { 
00145       cpl_msg_error(_id,"could not load arc image");
00146       /* destroy_image(im); */
00147       finddist_cfg_destroy (cfg);
00148       cpl_frameset_delete(stk);
00149       return -1;
00150     }
00151 
00152     /* load the fake on - fake off frame--- */
00153     imonind = load_image(cfg->nsFrame);
00154     if (imonind == NULL) { 
00155       cpl_msg_error(_id,"could not load on-off fake image");
00156       /* destroy_image(imonind); */
00157       destroy_image(im);
00158       finddist_cfg_destroy (cfg);
00159       cpl_frameset_delete(stk);
00160 
00161       return -1;
00162     }
00163 
00164     impoly  = load_image(cfg->nsFrame);
00165 
00166     if (impoly == NULL) { 
00167       cpl_msg_error(_id,"could not load on-off fake image");
00168       destroy_image(impoly);
00169       destroy_image(imonind);
00170       destroy_image(im);
00171       finddist_cfg_destroy (cfg);
00172       cpl_frameset_delete(stk);
00173       return -1;
00174     }
00175 
00176     /* load the fake on - fake off frame--- */
00177     mask = load_image(cfg->mask);
00178 
00179     if (mask == NULL) { 
00180       cpl_msg_error(_id,"could not load mask image");
00181       destroy_image(mask);
00182       destroy_image(impoly);
00183       destroy_image(imonind);
00184       destroy_image(im);
00185       finddist_cfg_destroy (cfg);
00186       cpl_frameset_delete(stk);
00187       return -1;
00188     }
00189 
00190 
00191     if( 0 != mul_img_local (im,mask) ) {
00192       cpl_msg_error(_id,"Failing to correct arc ima by bp mask");
00193       destroy_image(mask);
00194       destroy_image(impoly);
00195       destroy_image(imonind);
00196       destroy_image(im);
00197       finddist_cfg_destroy (cfg);
00198       cpl_frameset_delete(stk);
00199       return -1;
00200     }
00201 
00202 
00203     if( 0 != mul_img_local (imonind,mask) ) {
00204       cpl_msg_error(_id,"Failing to correct on-off fake ima by bp mask");
00205       destroy_image(mask);
00206       destroy_image(impoly);
00207       destroy_image(imonind);
00208       destroy_image(im);
00209       finddist_cfg_destroy (cfg);
00210       cpl_frameset_delete(stk);
00211       return -1;
00212     }
00213 
00214 
00215     if( 0 != mul_img_local (impoly,mask) ) {
00216       cpl_msg_error(_id,"Failing to correct on-off fake ima by bp mask");
00217       destroy_image(mask);
00218       destroy_image(impoly);
00219       destroy_image(imonind);
00220       destroy_image(im);
00221       finddist_cfg_destroy (cfg);
00222       cpl_frameset_delete(stk);
00223       return -1;
00224     }
00225     destroy_image(mask);
00226 
00227 
00228      /*
00229     lx = OneImage_lx_get( im );
00230     ly = OneImage_ly_get( im );
00231     */
00232 
00233     lx = im->lx;
00234     ly = im->ly;
00235 
00236 
00237     /* TO BE CHENGED THE FOLLOWING INPUT */
00238     /* open the line list and read the number of lines */
00239 
00240     tbl_line_list = cpl_table_load(cfg->lineList,1,0);
00241     if(cpl_error_get_code() != CPL_ERROR_NONE) {
00242       cpl_msg_error(_id,"problems loading table %s \n",cfg->lineList);
00243       cpl_msg_error(_id,(char* ) cpl_error_get_message());
00244  
00245       destroy_image(impoly);
00246       destroy_image(imonind);
00247       destroy_image(im);
00248       finddist_cfg_destroy (cfg);
00249       cpl_frameset_delete(stk);
00250 
00251      return -1;
00252     }
00253 
00254 
00255     n_lines = cpl_table_get_nrow(tbl_line_list);
00256 
00257  
00258     /*
00259     wave   = new_float_array (n_lines);
00260     intens = new_float_array (n_lines);
00261 
00262     if ( wave == NULL ||  intens == NULL ) {
00263       cpl_msg_error(_id, "could not allocate memory for the line list values" );
00264       cpl_table_delete(tbl_line_list);
00265       destroy_image(impoly);
00266       destroy_image(imonind);
00267       destroy_image(im);
00268       finddist_cfg_destroy (cfg);
00269       cpl_frameset_delete(stk);
00270       return -1;
00271     }
00272     */
00273 
00274 
00275 
00276 
00277 
00278     /* read the line list */
00279     /*
00280     n = cpl_table_get_nrow(tbl_line_list);
00281     n_lines = n;
00282     if(cpl_error_get_code() != CPL_ERROR_NONE) {
00283       cpl_msg_error(_id,(char* ) cpl_error_get_message());
00284       cpl_table_delete(tbl_line_list);
00285       destroy_image(impoly);
00286       destroy_image(imonind);
00287       destroy_image(im);
00288       finddist_cfg_destroy (cfg);
00289       cpl_frameset_delete(stk);
00290       return -1;
00291     }
00292     */
00293 
00294 
00295 
00296 
00297     /* THIS ORIGINATES A MEMORY LEAK 
00298     wave   = new_float_array (n);
00299     intens = new_float_array (n);
00300    
00301     if (wave == NULL || intens == NULL) {
00302       cpl_msg_error(_id,"could not allocate memory for the line list values\n" );
00303       return -1;
00304     }
00305     */
00306 
00307 
00308     wave   = cpl_table_get_data_float(tbl_line_list,"wave");
00309     if(cpl_error_get_code() != CPL_ERROR_NONE) {
00310       cpl_msg_error(_id,(char* ) cpl_error_get_message());
00311       cpl_table_delete(tbl_line_list);
00312       destroy_image(impoly);
00313       destroy_image(imonind);
00314       destroy_image(im);
00315       finddist_cfg_destroy (cfg);
00316       cpl_frameset_delete(stk);
00317 
00318       return -1;
00319     }
00320 
00321     intens = cpl_table_get_data_float(tbl_line_list,"int");
00322     if(cpl_error_get_code() != CPL_ERROR_NONE) {
00323       cpl_msg_error(_id,(char* ) cpl_error_get_message());
00324       cpl_table_delete(tbl_line_list);
00325       destroy_image(impoly);
00326       destroy_image(imonind);
00327       destroy_image(im);
00328       finddist_cfg_destroy (cfg);
00329       cpl_frameset_delete(stk);
00330       return -1;
00331     }
00332 
00333 
00334 
00335     /*
00336       destroy_array ( wave );
00337       destroy_array ( intens );
00338     */
00339 
00340     /* 
00341        #---------------------------------------------------------------------
00342        #---------------------------FINDLINES---------------------------------
00343        #---------------------------------------------------------------------
00344     */
00345 
00346     cpl_msg_info(_id,"Find Lines");
00347     /* do the wavelength calibration, that means: 
00348        find the dispersion relation and parameterize its coefficients
00349     */
00350 
00351     acoefs  = new_2Dfloat_array(cfg->nrDispCoefficients, lx);
00352     /* allocate memory */
00353     n_found_lines    = new_int_array(lx); 
00354     row_clean        = new_2Dint_array(lx, n_lines);
00355     wavelength_clean = new_2Dfloat_array(lx, n_lines);
00356     sum_pointer      = new_int_array(1);
00357  
00358     /* find the emission lines in each image column */
00359     intarray_set_value(sum_pointer, 0, 0);
00360 
00361     check = findLines(im, 
00362                       wave, 
00363                       intens, 
00364                       n_lines, 
00365                       row_clean, 
00366                       wavelength_clean, 
00367                       cfg->guessBeginWavelength, 
00368                       cfg->guessDispersion1,
00369                       cfg->guessDispersion2, 
00370                       cfg->mindiff, 
00371                       cfg->halfWidth, 
00372                       n_found_lines, 
00373                       cfg->sigma, 
00374                       sum_pointer );
00375 
00376 
00377 
00378    cpl_table_delete(tbl_line_list);
00379 
00380     if (-1 == check ) { 
00381       cpl_msg_error (_id,"FindLines failed!");
00382     /* deallocate memory */
00383  
00384     destroy_intarray ( n_found_lines );
00385     destroy_2Dintarray (row_clean, lx);
00386     destroy_2Darray ( wavelength_clean, lx );
00387     destroy_intarray ( sum_pointer );
00388     destroy_2Darray ( acoefs, cfg->nrDispCoefficients );
00389 
00390 
00391       destroy_image(impoly);
00392       destroy_image(imonind);
00393       destroy_image(im);
00394       finddist_cfg_destroy (cfg);
00395       cpl_frameset_delete(stk);
00396 
00397       return -1;
00398     }
00399 
00400  
00401     /*
00402 #----------------------------------------------------------------------------
00403 #------------------------------WAVECALIB-------------------------------------
00404 #----------------------------------------------------------------------------
00405     */
00406     cpl_msg_info(_id,"Do wave calib");
00407     sum = intarray_get_value(sum_pointer,0);
00408 
00409     /* #allocate memory for the fit parameters */
00410     par = newFitParams( sum );
00411     if ( par == NULL ) {
00412       cpl_msg_error(_id,"newFitParams failed!");
00413       /* destroyFitParams(par); */
00414       destroy_intarray ( n_found_lines );
00415       destroy_2Dintarray (row_clean, lx);
00416       destroy_2Darray ( wavelength_clean, lx );
00417       destroy_intarray ( sum_pointer );
00418       destroy_2Darray ( acoefs, cfg->nrDispCoefficients );
00419       destroy_image(im);
00420       destroy_image(imonind);
00421       destroy_image(impoly);
00422       finddist_cfg_destroy (cfg);
00423       cpl_frameset_delete(stk);
00424 
00425       return -1;
00426     }
00427 
00428 
00429     /* fit each line, make a polynomial fit and fit the resulting fit 
00430        coefficients across the columns of the slitlet
00431     */
00432 map = waveCal(im, 
00433               par, 
00434               acoefs, 
00435               cfg->nslitlets, 
00436               row_clean, 
00437               wavelength_clean,
00438               n_found_lines, 
00439               cfg->guessDispersion1, 
00440               cfg->halfWidth, 
00441               cfg->minAmplitude, 
00442               cfg->maxResidual, 
00443               cfg->fwhm, 
00444               cfg->nrDispCoefficients, 
00445               cfg->nrCoefCoefficients, 
00446               cfg->sigmaFactor, 
00447               cfg->pixeldist, 
00448               cfg->pixel_tolerance);
00449 
00450 
00451  if (NULL == map ) { 
00452    cpl_msg_error (_id,"waveCal failed!");
00453    destroyFitParams(par); 
00454    destroy_intarray ( n_found_lines );
00455    destroy_2Dintarray (row_clean, lx);
00456    destroy_2Darray ( wavelength_clean, lx );
00457    destroy_intarray ( sum_pointer );
00458    destroy_2Darray ( acoefs, cfg->nrDispCoefficients );
00459    destroy_image(impoly);
00460    destroy_image(imonind);
00461    destroy_image(im);
00462    cpl_frameset_delete(stk);
00463    finddist_cfg_destroy (cfg);
00464    return -1;
00465  }
00466 
00467 
00468 
00469  shift = checkLinePositions (im, acoefs, cfg->nrDispCoefficients, par);
00470 
00471  if (FLAG == shift) { 
00472    cpl_msg_error (_id, "checkForLinePositions failed!");
00473    destroy_2Darray ( wavelength_clean, lx );
00474    destroy_2Dintarray (row_clean, lx);
00475    destroy_intarray ( n_found_lines );
00476    destroy_intarray ( sum_pointer );
00477    destroyFitParams(par);
00478    destroy_image(impoly);
00479    destroy_image(imonind);
00480    destroy_image(im);
00481    cpl_frameset_delete(stk);
00482    finddist_cfg_destroy (cfg);
00483    return -1;
00484  }
00485  destroy_image(map);
00486 
00487 
00488 
00489 
00490  /* free memory */
00491  destroy_2Darray ( wavelength_clean, lx );
00492  destroy_2Dintarray (row_clean, lx);
00493  destroy_intarray ( n_found_lines );
00494  destroy_intarray ( sum_pointer );
00495  destroy_2Darray ( acoefs, cfg->nrDispCoefficients );
00496   /* To not have warning
00497  destroy_array ( wave );
00498  destroy_array ( intens );
00499  */
00500 
00501 
00502 
00503  /*
00504 #----------------------------------------------------------------------------
00505 #-------------------------SLITFITS-------------------------------------------
00506 #----------------------------------------------------------------------------
00507 #--fit the slitlet edge positions if desired--
00508  */
00509 
00510 
00511 
00512 
00513  cpl_msg_info(_id,"Find slit pos");
00514  /*allocate memory for the slitlet position array */
00515  slit_pos = new_2Dfloat_array(32,2);
00516 
00517  fit = fitSlitsBoltz( im, par, slit_pos, cfg->boxLength, 
00518                       cfg->yBox, cfg->diffTol );
00519 
00520  if (fit < 0) {
00521    cpl_msg_error(_id, "fitSlitsBoltz failed" );
00522    destroy_2Darray ( slit_pos, 32 );
00523    destroy_image(im);
00524    destroy_image(impoly);
00525    destroy_image(imonind);
00526    destroyFitParams(par); 
00527    return -1;
00528  }
00529  destroyFitParams(par); 
00530 
00531  tbl_spos=cpl_table_new(32);
00532  cpl_table_new_column(tbl_spos,"pos1", CPL_TYPE_DOUBLE);
00533  cpl_table_new_column(tbl_spos,"pos2", CPL_TYPE_DOUBLE);
00534  cpl_table_set_column_format(tbl_spos,"pos1", "15.9f");
00535  cpl_table_set_column_format(tbl_spos,"pos2", "15.9f");
00536 
00537  for (i =0; i< 32; i++) {
00538      /*
00539      fprintf( slitpos_file, "%15.9f %15.9f \n",
00540       array2D_get_value(slit_pos,i,0),array2D_get_value(slit_pos,i,1));
00541      */
00542       cpl_table_set_double(tbl_spos,"pos1",i,array2D_get_value(slit_pos,i,0));
00543       cpl_table_set_double(tbl_spos,"pos2",i,array2D_get_value(slit_pos,i,1));
00544      
00545  }
00546  strcpy(tbl_name,"out_slitlets_pos_predist.tfits");
00547    
00548 
00549     if(-1 == sinfoni_pro_dump_tbl(tbl_spos,stk,sof,tbl_name,
00550          PRO_SLITLETS_POS_PREDIST,plugin_id,config)) {
00551          cpl_msg_error(_id,"cannot dump tbl %s", tbl_name);
00552      destroy_2Darray ( slit_pos, 32 );
00553      destroy_image(im);
00554      destroy_image(impoly);
00555      destroy_image(imonind);
00556          cpl_table_delete(tbl_spos);
00557          cpl_frameset_delete(stk);
00558          finddist_cfg_destroy (cfg);
00559          return -1;
00560     }
00561    
00562     cpl_table_delete(tbl_spos);
00563 
00564  /*
00565 #---------------------------------------------------------
00566 # do the north - south - test
00567 #---------------------------------------------------------
00568 #distances = north_south_test(imonind, nslits, nshalfWidth, nsfwhm , minDiff, estimatedDist, devtol,512,1536 )
00569  */
00570  cpl_msg_info(_id,"Do north south test");
00571 
00572  cpl_msg_info(_id,"cfg->nslits =%d\n", cfg->nslits);
00573  cpl_msg_info(_id,"cfg->nshalfWidth =%d\n", cfg->nshalfWidth);
00574  cpl_msg_info(_id,"cfg->nsfwhm=%f\n",cfg->nsfwhm ); 
00575  cpl_msg_info(_id,"cfg->minDiff=%f\n", cfg->minDiff);
00576  cpl_msg_info(_id,"cfg->estimated_dist=%f\n", cfg->estimated_dist);
00577  cpl_msg_info(_id,"cfg->devtol=%f\n", cfg->devtol);
00578  cpl_msg_info(_id,"cfg->loPos=%d\n", cfg->loPos);
00579  cpl_msg_info(_id,"cfg->hiPos=%d\n",cfg->hiPos);
00580 
00581 distances = north_south_test(imonind, 
00582                              cfg->nslits, 
00583                              cfg->nshalfWidth, 
00584                              cfg->nsfwhm , 
00585                              cfg->minDiff, 
00586                              cfg->estimated_dist, 
00587                              cfg->devtol,
00588                              cfg->loPos,
00589                              cfg->hiPos);
00590 
00591 
00592       destroy_image(imonind); 
00593 
00594 
00595 
00596  if (distances == NULL ){
00597    cpl_msg_error (_id,"north_south_test failed");
00598  
00599    /*destroy_array(distances); */
00600      destroy_2Darray ( slit_pos, 32 );
00601       destroy_image(im);
00602       destroy_image(impoly);
00603       cpl_frameset_delete(stk);
00604       finddist_cfg_destroy (cfg);
00605 
00606   return -1;
00607  }
00608 
00609 
00610     
00611  /*
00612 #---------------------------------------------------------
00613 # get firstcol
00614 #---------------------------------------------------------
00615  */
00616  cpl_msg_info(_id,"get first col");
00617  first = new_float_array(61);
00618  total_dist=0;
00619  n=0;
00620 
00621 
00622  /*
00623  cpl_msg_info(_id,"%f %f %f\n", 
00624        array2D_get_value(slit_pos,0,0), 
00625        array2D_get_value(slit_pos,0,1), 
00626        array2D_get_value(slit_pos,0,1)-array2D_get_value(slit_pos,0,0);
00627  */
00628 
00629  for (i=0; i< 31; i++) {
00630    total_dist=total_dist+array_get_value(distances,i);
00631 
00632    /*
00633     cpl_msg_info(_id,"%f %f %f\n", 
00634             array2D_get_value(slit_pos,i+1,0), 
00635             array2D_get_value(slit_pos,i+1,1), 
00636             array2D_get_value(slit_pos,i+1,1)-
00637             array2D_get_value(slit_pos,i+1,0));
00638    */
00639    
00640 
00641 
00642    for (j=0; j< 2; j++) {
00643      if (j == 0) {
00644        if (i != 30) {
00645      newval = array2D_get_value(slit_pos,i+1,j) - total_dist;
00646      array_set_value(first, newval, n);
00647      n = n+1;
00648        }
00649      }
00650      if (j == 1) {
00651        newval = array2D_get_value(slit_pos,i,j) - total_dist;
00652        array_set_value(first, newval, n);
00653        n = n+1;
00654      }
00655    }
00656  }
00657  value = f_median(first,61);
00658  cpl_msg_info(_id,"Firstcol =%f", value);
00659  destroy_array(first);
00660  destroy_array(distances);
00661  destroy_2Darray ( slit_pos, 32 );
00662 
00663 
00664 
00665  
00666    /*
00667 #---------------------------------------------------------
00668 # Determine distortions
00669 #---------------------------------------------------------
00670    */
00671  cpl_msg_info(_id,"Determine distortions");
00672  degx=cpl_calloc(8,sizeof(int));
00673  degy=cpl_calloc(8,sizeof(int));
00674  coef=cpl_calloc(8,sizeof(double));
00675 
00676 
00677 
00678  /* this generates a mem leak */
00679 
00680  poly=compute_distortion(impoly,0,lx,0,ly,1,value);
00681 
00682  /* up to here ok */
00683  /* To check mem leaks */
00684  destroy_image(impoly);
00685 
00686 
00687  /*
00688  if(poly==NULL) {
00689    cpl_msg_error(_id,"Error computing distortions");
00690    destroy_image(impoly);
00691    cpl_free(degx);
00692    cpl_free(degy);
00693    cpl_free(coef);
00694    cpl_free(poly_file);
00695    destroy_image(im);
00696    cpl_frameset_delete(stk);
00697    finddist_cfg_destroy (cfg);
00698    return -1;
00699  }
00700  */
00701 
00702  cpl_msg_info(_id,"Save results");
00703  poly_file = (char*) cpl_calloc(MAX_NAME_SIZE,sizeof(char*));
00704 
00705  strcpy(poly_file,"tmp_polt.list"); 
00706  poly2d_save(poly,poly_file);
00707 
00708 
00709 
00710 if((poly_in=fopen(poly_file,"r"))==NULL) {
00711         cpl_msg_error("Cannot open file %s", poly_file);
00712         return -1;
00713 }
00714 fscanf(poly_in,"%s %s %s %s",temp,temp,temp,temp);
00715 ch=fgetc(poly_in);
00716 while(ch!='p'){
00717     ch=fgetc(poly_in);
00718     poly_string[count]=ch;
00719     count++;
00720 }
00721 
00722 
00723 
00724 
00725 poly_string[count-2]='\0';
00726 /*fscanf(poly_in,"%s",poly_string);*/
00727 cpl_msg_warning("%s",poly_string);
00728 
00729   sscanf(poly_string,"%d %d %f" 
00730                      "%d %d %f" 
00731                      "%d %d %f" 
00732                      "%d %d %f" 
00733                      "%d %d %f"
00734                      "%d %d %f" 
00735                      "%d %d %f" 
00736                      "%d %d %f",
00737                  &dx0,&dy0,&c00,
00738                  &dx1,&dy0,&c10,
00739                  &dx0,&dy1,&c01,
00740                  &dx1,&dy1,&c11,
00741                  &dx2,&dy1,&c20,
00742                  &dx1,&dy2,&c02,
00743                  &dx2,&dy2,&c21,
00744                  &dx1,&dy2,&c12);
00745 
00746 
00747 
00748   poly_tbl=cpl_table_new(8);
00749   cpl_table_new_column(poly_tbl,"degx", CPL_TYPE_INT);
00750   cpl_table_new_column(poly_tbl,"degy", CPL_TYPE_INT);
00751   cpl_table_new_column(poly_tbl,"coeff", CPL_TYPE_DOUBLE);
00752   cpl_table_set_int(poly_tbl,"degx",0,0);
00753   cpl_table_set_int(poly_tbl,"degx",1,1);
00754   cpl_table_set_int(poly_tbl,"degx",2,0);
00755   cpl_table_set_int(poly_tbl,"degx",3,1);
00756   cpl_table_set_int(poly_tbl,"degx",4,2);
00757   cpl_table_set_int(poly_tbl,"degx",5,0);
00758   cpl_table_set_int(poly_tbl,"degx",6,2);
00759   cpl_table_set_int(poly_tbl,"degx",7,1);
00760 
00761   cpl_table_set_int(poly_tbl,"degy",0,0);
00762   cpl_table_set_int(poly_tbl,"degy",1,0);
00763   cpl_table_set_int(poly_tbl,"degy",2,1);
00764   cpl_table_set_int(poly_tbl,"degy",3,1);
00765   cpl_table_set_int(poly_tbl,"degy",4,0);
00766   cpl_table_set_int(poly_tbl,"degy",5,2);
00767   cpl_table_set_int(poly_tbl,"degy",6,1);
00768   cpl_table_set_int(poly_tbl,"degy",7,2);
00769 
00770 
00771   cpl_table_set_double(poly_tbl,"coeff",0,c00);
00772   cpl_table_set_double(poly_tbl,"coeff",1,c10);
00773   cpl_table_set_double(poly_tbl,"coeff",2,c01);
00774   cpl_table_set_double(poly_tbl,"coeff",3,c11);
00775   cpl_table_set_double(poly_tbl,"coeff",4,c20);
00776   cpl_table_set_double(poly_tbl,"coeff",5,c02);
00777   cpl_table_set_double(poly_tbl,"coeff",6,c21);
00778   cpl_table_set_double(poly_tbl,"coeff",7,c12);
00779 
00780   cpl_msg_info(_id,"Polynomial distortion coefficients \n"); 
00781   cpl_msg_info(_id,"c= %g   %g  %g  %g  %g  %g  %g  %g\n",
00782                    c00,c10,c01,c11,c20,c02,c21,c12);
00783 
00784 
00785 
00786 
00787 
00788        /* QC LOG */
00789 
00790        
00791   qclog_tbl = sinfoni_qclog_init(13);
00792   key_name=cpl_calloc(FILE_NAME_SZ,sizeof(char));
00793   key_value=cpl_calloc(FILE_NAME_SZ,sizeof(char));
00794 
00795   sprintf(key_name,"%s%d%d","QC COEFF",0,0);
00796   sprintf(key_value,"%g",c00);
00797   sinfoni_qclog_add(qclog_tbl,0,key_name,"CPL_TYPE_DOUBLE",
00798             key_value,"Polynomial distortion coefficient");
00799 
00800   sprintf(key_name,"%s%d%d","QC COEFF",1,0);
00801   sprintf(key_value,"%g",c10);
00802   sinfoni_qclog_add(qclog_tbl,1,key_name,"CPL_TYPE_DOUBLE",
00803             key_value,"Polynomial distortion coefficient");
00804 
00805   sprintf(key_name,"%s%d%d","QC COEFF",0,1);
00806   sprintf(key_value,"%g",c01);
00807   sinfoni_qclog_add(qclog_tbl,2,key_name,"CPL_TYPE_DOUBLE",
00808             key_value,"Polynomial distortion coefficient");
00809 
00810   sprintf(key_name,"%s%d%d","QC COEFF",1,1);
00811   sprintf(key_value,"%g",c11);
00812   sinfoni_qclog_add(qclog_tbl,3,key_name,"CPL_TYPE_DOUBLE",
00813             key_value,"Polynomial distortion coefficient");
00814 
00815   sprintf(key_name,"%s%d%d","QC COEFF",2,0);
00816   sprintf(key_value,"%g",c20);
00817   sinfoni_qclog_add(qclog_tbl,4,key_name,"CPL_TYPE_DOUBLE",
00818             key_value,"Polynomial distortion coefficient");
00819 
00820   sprintf(key_name,"%s%d%d","QC COEFF",0,2);
00821   sprintf(key_value,"%g",c02);
00822   sinfoni_qclog_add(qclog_tbl,5,key_name,"CPL_TYPE_DOUBLE",
00823             key_value,"Polynomial distortion coefficient");
00824 
00825   sprintf(key_name,"%s%d%d","QC COEFF",2,1);
00826   sprintf(key_value,"%g",c21);
00827   sinfoni_qclog_add(qclog_tbl,6,key_name,"CPL_TYPE_DOUBLE",
00828             key_value,"Polynomial distortion coefficient");
00829 
00830   sprintf(key_name,"%s%d%d","QC COEFF",1,2);
00831   sprintf(key_value,"%g",c12);
00832   sinfoni_qclog_add(qclog_tbl,7,key_name,"CPL_TYPE_DOUBLE",
00833             key_value,"Polynomial distortion coefficient");
00834 
00835 
00836 
00837   sprintf(key_name,"%s","QC OFFSET");
00838   sprintf(key_value,"%g",value);
00839   sinfoni_qclog_add(qclog_tbl,8,key_name,"CPL_TYPE_DOUBLE",
00840             key_value,"Polynomial distortion offset");
00841 
00842 
00843   /*
00844   for(i=0;i<8;i++) {
00845 
00846     printf("tc = %g\n", cpl_table_get_double(poly_tbl,"coeff",i,status));
00847 
00848   } 
00849  
00850   find_distortion_qclog(poly_tbl, qclog_tbl);
00851   */
00852  
00853   xshift=compute_shift(x_c,y_c,c00,c10,c01,c11,c20,c02,c21,c12);
00854   sprintf(key_value,"%g",xshift);
00855   sinfoni_qclog_add(qclog_tbl,9,"QC XSHIFT CC","CPL_TYPE_DOUBLE",
00856             key_value,"X shift in x_c,y_c");
00857 
00858 
00859   xshift=compute_shift(x_l,y_l,c00,c10,c01,c11,c20,c02,c21,c12);
00860   sprintf(key_value,"%g",xshift);
00861   sinfoni_qclog_add(qclog_tbl,10,"QC XSHIFT LL","CPL_TYPE_DOUBLE",
00862             key_value,"X shift in x_l,y_l");
00863 
00864 
00865   xshift=compute_shift(x_l,y_u,c00,c10,c01,c11,c20,c02,c21,c12);
00866   sprintf(key_value,"%g",xshift);
00867   sinfoni_qclog_add(qclog_tbl,8,"QC XSHIFT UL","CPL_TYPE_DOUBLE",
00868             key_value,"X shift in x_l,y_u");
00869 
00870   xshift=compute_shift(x_u,y_u,c00,c10,c01,c11,c20,c02,c21,c12);
00871   sprintf(key_value,"%g",xshift);
00872   sinfoni_qclog_add(qclog_tbl,11,"QC XSHIFT UR","CPL_TYPE_DOUBLE",
00873             key_value,"X shift in x_u,y_u");
00874 
00875 
00876   xshift=compute_shift(x_u,y_l,c00,c10,c01,c11,c20,c02,c21,c12);
00877   sprintf(key_value,"%g",xshift);
00878   sinfoni_qclog_add(qclog_tbl,12,"QC XSHIFT LR","CPL_TYPE_DOUBLE",
00879             key_value,"X shift in x_u,y_l");
00880 
00881 
00882  
00883  
00884   if(-1 == sinfoni_pro_save_tbl(poly_tbl,stk,sof,cfg->outName,
00885          PRO_DISTORTION,qclog_tbl,plugin_id,config)) {
00886          cpl_msg_error(_id,"cannot dump tbl %s", cfg->outName);
00887      cpl_table_delete(poly_tbl);
00888      cpl_table_delete(qclog_tbl);
00889 
00890 
00891      cpl_free(key_name);
00892      cpl_free(key_value);
00893      cpl_free(degx);
00894      cpl_free(degy);
00895      cpl_free(coef);
00896      cpl_free(poly_file);
00897      destroy_image(im);
00898      poly2d_free(poly);
00899      cpl_frameset_delete(stk);
00900      finddist_cfg_destroy (cfg);
00901 
00902 
00903          return -1;
00904   }
00905   cpl_table_delete(poly_tbl);
00906   cpl_table_delete(qclog_tbl);
00907 
00908 
00909   cpl_free(key_name);
00910   cpl_free(key_value);
00911   cpl_free(degx);
00912   cpl_free(degy);
00913   cpl_free(coef);
00914   cpl_free(poly_file);
00915   destroy_image(im);
00916   poly2d_free(poly);
00917   cpl_frameset_delete(stk);
00918   finddist_cfg_destroy (cfg);
00919 
00920     return 0;
00921 
00922 }
00923 
00924 
00925 
00926 int
00927 find_distortion_qclog(cpl_table* poly_tbl, cpl_table* qclog_tbl)
00928 {
00929     char* key_value=NULL;
00930     char* key_name=NULL;
00931     int* status=NULL;
00932 
00933     int i=0;
00934     int n=0;
00935     int degx=0;
00936     int degy=0;
00937     double coef=0;
00938 
00939 
00940     /* QC LOG */
00941     n=cpl_table_get_nrow(poly_tbl);
00942     qclog_tbl = sinfoni_qclog_init(n);
00943 
00944     key_name  = cpl_calloc(MAX_NAME_SIZE,sizeof(char));
00945     key_value = cpl_calloc(MAX_NAME_SIZE,sizeof(char));
00946     for (i=0;i<n;i++) {
00947         degx=cpl_table_get_int(poly_tbl,"degx",i,status);
00948         degy=cpl_table_get_int(poly_tbl,"degy",i,status);
00949         coef=cpl_table_get_double(poly_tbl,"coeff",i,status);
00950         sprintf(key_name,"%s%d%d","QC COEFF",degx,degy);
00951         sprintf(key_value,"%g",coef);
00952         sinfoni_qclog_add(qclog_tbl,i,key_name,"CPL_TYPE_DOUBLE",
00953               key_value,"Polynomial distortion coefficient");
00954 
00955     }
00956 
00957     cpl_free(key_name);
00958     cpl_free(key_value);
00959     return 0;
00960 
00961 }
00962 
00963 double
00964 compute_shift(double x,
00965               double y,
00966               double c00,
00967               double c10,
00968               double c01,
00969               double c11,
00970               double c20,
00971               double c02,
00972               double c21,
00973               double c12)
00974 {
00975 
00976   double x_shift=0;
00977   double shift=0;
00978   
00979   shift=c00+c10*x+c01*y+c11*x*y+c20*x*x+c02*y*y+c21*x*x*y+c12*x*y*y;
00980   x_shift=x-shift;
00981   return x_shift;
00982 
00983 
00984 
00985 }

Generated on Wed Oct 26 13:08:52 2005 for SINFONI Pipeline Reference Manual by doxygen1.2.13.1 written by Dimitri van Heesch, © 1997-2001