sinfo_bp_lin.c

00001 /*
00002  * This file is part of the ESO SINFONI Pipeline
00003  * Copyright (C) 2004,2005 European Southern Observatory
00004  *
00005  * This program is free software; you can redistribute it and/or modify
00006  * it under the terms of the GNU General Public License as published by
00007  * the Free Software Foundation; either version 2 of the License, or
00008  * (at your option) any later version.
00009  *
00010  * This program is distributed in the hope that it will be useful,
00011  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00012  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00013  * GNU General Public License for more details.
00014  *
00015  * You should have received a copy of the GNU General Public License
00016  * along with this program; if not, write to the Free Software
00017  * Foundation, 51 Franklin St, Fifth Floor, Boston, MA  02111-1307  USA
00018  */
00019 /*----------------------------------------------------------------------------
00020    
00021    File name    :       sinfo_bp_lin.c
00022    Author       :   J. Schreiber
00023    Created on   :   May 5, 2003
00024    Description  :   Different methods for searching for bad pixels
00025                         used in the recipe spiffi_badsearch 
00026 
00027  ---------------------------------------------------------------------------*/
00028 #ifdef HAVE_CONFIG_H
00029 #  include <config.h>
00030 #endif
00031 
00032 /*----------------------------------------------------------------------------
00033                                 Includes
00034  ---------------------------------------------------------------------------*/
00035 //#include <cpl_imagelist_basic.h>
00036 #include <cpl.h>
00037 #include "irplib_calib.h"
00038 #include "sinfo_bp_lin.h"
00039 #include "sinfo_detlin.h"
00040 #include "sinfo_detlin_ini_by_cpl.h"
00041 #include "sinfo_pro_save.h"
00042 #include "sinfo_pro_types.h"
00043 #include "sinfo_functions.h"
00044 #include "sinfo_hidden.h"
00045 #include "sinfo_error.h"
00046 #include "sinfo_utils_wrappers.h"
00047 
00048 
00055 /*----------------------------------------------------------------------------
00056                                 Defines
00057  ---------------------------------------------------------------------------*/
00058 /*----------------------------------------------------------------------------
00059                              Function Definitions
00060  ---------------------------------------------------------------------------*/
00087 int 
00088 sinfo_new_bp_search_lin(const char* plugin_id, 
00089                         cpl_parameterlist* config, 
00090                         cpl_frameset* sof)
00091 {
00092     detlin_config * cfg=NULL;
00093     cpl_imagelist * imgl=NULL ;
00094     cpl_imagelist * coeffs_imgl=NULL;
00095     cpl_image * img_tmp=NULL ;
00096     cpl_image * mask=NULL;
00097     cpl_vector* cube_mean=NULL;
00098     cpl_frameset* raw=NULL;
00099     cpl_frameset* raw_on=NULL;
00100     cpl_frameset* raw_of=NULL;
00101     cpl_frame* frm=NULL;
00102     cpl_frame* frm_dup=NULL;
00103 
00104     cpl_table* qclog_tbl=NULL;
00105 
00106     cpl_parameter *p=NULL;             
00107     cpl_polynomial *pol=NULL;
00108     cpl_vector* vec_adl=NULL;
00109     cpl_vector* vec_med=NULL;
00110     double* mse=NULL;
00111     int i=0;
00112     int n_bad=0 ;
00113     int nff=0;
00114     char key_value[FILE_NAME_SZ];
00115     char key_name[FILE_NAME_SZ];
00116     cpl_table* det_lin=NULL;
00117     int nraw=0;
00118     int* status=NULL;
00119     cpl_table* gain=NULL;
00120     int ngain=0;
00121     int sz_imgl=0;
00122 
00123      int llx=270;
00124      int lly=1000;
00125      int urx=320;
00126      int ury=1050;
00127      int zone1[4];
00128      int zone2[4];
00129      //int kappa=3;
00130      //int nclip=5;
00131 
00132 
00133     /* parse the file names and parameters to the bad_config 
00134        data structure cfg */
00135     check_nomsg(raw=cpl_frameset_new());
00136     cknull(cfg = sinfo_parse_cpl_input_detlin(config,sof,&raw),
00137            "could not parse .ini file!");
00138 
00139     /* ======================================================================
00140        DETERMINES LINEARITY COEFF AS DFO DOES
00141        ======================================================================
00142     */
00143     nff=cpl_frameset_get_size(raw);
00144     raw_on=cpl_frameset_new();
00145     raw_of=cpl_frameset_new();
00146     /* separates on and off frames */
00147 
00148     for(i=0;i<nff;i++) {
00149       frm=cpl_frameset_get_frame(raw,i);
00150       frm_dup=cpl_frame_duplicate(frm);
00151       if(sinfo_frame_is_on(frm)) {
00152     cpl_frameset_insert(raw_on,frm_dup);
00153       } else {
00154     cpl_frameset_insert(raw_of,frm_dup);
00155       }
00156     }  
00157     /*
00158     cknull(det_lin=sinfo_get_linearity(raw_on,raw_of),
00159            "Error computing linearity");
00160     */
00161     zone1[0]=llx;
00162     zone1[1]=lly;
00163     zone1[2]=urx;
00164     zone1[3]=ury;
00165     zone2[0]=270;
00166     zone2[1]=1030;
00167     zone2[2]=310;
00168     zone2[3]=1060;
00169 
00170     /*
00171     cknull(det_lin=irplib_compute_linearity(raw_on,raw_of),
00172            "Error computing linearity");
00173     */
00174     sinfo_msg("computes linearity");
00175     cknull(det_lin=sinfo_compute_linearity(raw_on,raw_of),
00176            "Error computing linearity");
00177 
00178     check_nomsg(nraw=cpl_table_get_nrow(det_lin));
00179     check_nomsg(vec_adl=cpl_vector_new(nraw));
00180     check_nomsg(vec_med=cpl_vector_new(nraw));
00181 
00182     for(i=0;i<nraw;i++) {
00183       cpl_vector_set(vec_adl,i,cpl_table_get_double(det_lin,"adl",i,status)); 
00184       cpl_vector_set(vec_med,i,cpl_table_get_double(det_lin,"med",i,status));
00185     }
00186     check_nomsg(pol=cpl_polynomial_fit_1d_create(vec_adl,vec_med,
00187                                                  cfg->order,mse));
00188     sinfoni_free_vector(&vec_adl);
00189     sinfoni_free_vector(&vec_med);
00190 
00191     cknull_nomsg(qclog_tbl=sinfo_qclog_init());
00192 
00193     for(i=0;i<cfg->order+1;i++) {
00194 
00195       snprintf(key_name,MAX_NAME_SIZE-1,"%s%i%s","QC BP-MAP LIN",i," MED");
00196       sinfo_qclog_add_double(qclog_tbl,key_name,
00197                              cpl_polynomial_get_coeff(pol,&i),
00198                              "Linearity Polynomial Coeff","%g");
00199     }
00200     sinfo_free_polynomial(&pol);
00201    
00202     check_nomsg(p = cpl_parameterlist_find(config, "sinfoni.bp.method"));
00203     snprintf(key_value, MAX_NAME_SIZE-1,"%s",cpl_parameter_get_string(p));
00204     ck0_nomsg(sinfo_qclog_add_string(qclog_tbl,"QC BP-MAP METHOD",key_value,
00205                 "BP search method","%s"));
00206 
00207     ck0(sinfo_pro_save_tbl(det_lin,raw,sof,BP_LIN_LIN_DET_INFO_OUT_FILENAME,
00208         "LIN_DET_INFO",qclog_tbl,plugin_id,config),
00209         "cannot dump ims %s",BP_LIN_LIN_DET_INFO_OUT_FILENAME);
00210 
00211     sinfo_free_table(&det_lin);
00212     sinfo_free_table(&qclog_tbl);
00213 
00214     /*=======================================================*/
00215     sinfo_msg("computes gain");
00216 
00217     cknull(gain=sinfo_compute_gain(raw_on,raw_of),"Error computing gain");
00218     /*
00219     cknull(gain=irplib_compute_gain(raw_on,raw_of,zone1,zone2,kappa,nclip),
00220              "Error computing gain");
00221     */
00222     sinfo_free_frameset(&raw_on);
00223     sinfo_free_frameset(&raw_of);
00224 
00225     check_nomsg(ngain=cpl_table_get_nrow(gain)); 
00226     cknull_nomsg(qclog_tbl=sinfo_qclog_init());
00227 
00228     for(i=0;i<ngain;i++) {
00229         snprintf(key_name,MAX_NAME_SIZE-1,"%s%i","QC GAIN",i);
00230         ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,
00231                                 cpl_table_get_double(gain,"gain",i,status),
00232                                 "Detector gain","%g"));
00233     }
00234 
00235     ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,"QC GAIN",
00236      cpl_table_get_column_median(gain,"gain"),"Detector gain median","%g"));
00237 
00238     ck0(sinfo_pro_save_tbl(gain,raw,sof,BP_LIN_GAIN_OUT_FILENAME,"GAIN_INFO",
00239         qclog_tbl,plugin_id,config),
00240         "cannot dump tbl %s", BP_LIN_GAIN_OUT_FILENAME);   
00241 
00242     sinfo_free_table(&gain);
00243     sinfo_free_table(&qclog_tbl);
00244 
00245     /* =======================================================================
00246        DETERMINES LINEARITY COEFF AS MPE DOES
00247        =======================================================================
00248     */
00249 
00250 
00251     check_nomsg(imgl = cpl_imagelist_new());
00252     for ( i = 0 ; i < cfg->nframes ; i++ )
00253     {
00254       if(sinfo_is_fits_file (cfg->framelist[i]) != 1) {
00255     sinfo_msg_error("Input file %s is not FITS",cfg->framelist[i] );
00256     return -1;
00257       }
00258       check_nomsg(img_tmp=cpl_image_load(cfg->framelist[i],CPL_TYPE_FLOAT,0,0));
00259       check_nomsg(cpl_imagelist_set(imgl,img_tmp,i));
00260     }
00261 
00262     /*----------------------------------------------------------------
00263      *---------------------- SEARCH FOR BAD PIXELS---------------------
00264      *--------------------------------------------------------------*/
00265     sinfo_msg("Search for bad pixels");
00266     cknull(coeffs_imgl=sinfo_new_fit_intensity_course(imgl, 
00267                                                       cfg->order, 
00268                               cfg->loReject, 
00269                                                       cfg->hiReject),
00270        "could not fit polynomial and store coeffs in a data cube!");
00271 
00272      sinfo_free_imagelist(&imgl) ;
00273     check_nomsg(sz_imgl=cpl_imagelist_get_size(coeffs_imgl));
00274 
00275 
00276     /*---store the polynomial fit coefficients in a data cube----*/
00277     /* CUBE NOT DUMPED BECAUSE PROBLEMS ON LAST PLANE */
00278     cube_mean=cpl_vector_new(sz_imgl);
00279 
00280      /* QC LOG */
00281     cknull_nomsg(qclog_tbl=sinfo_qclog_init());
00282     for(i=0;i<sz_imgl;i++) {
00283        cpl_vector_set(cube_mean,i,
00284                       cpl_image_get_median(cpl_imagelist_get(coeffs_imgl,i)));
00285 
00286        snprintf(key_name,MAX_NAME_SIZE-1,"%s%i%s","QC BP-MAP LIN",i," MEAN");
00287        snprintf(key_value,MAX_NAME_SIZE-1,"%g",cpl_vector_get(cube_mean,i));
00288        ck0_nomsg(sinfo_qclog_add_double(qclog_tbl,key_name,
00289                                         cpl_vector_get(cube_mean,i),
00290                                         "Linearity Polynomial Coeff","%g"));
00291 
00292     }
00293     sinfoni_free_vector(&cube_mean);
00294     check_nomsg(p = cpl_parameterlist_find(config, "sinfoni.bp.method"));
00295     snprintf(key_value, MAX_NAME_SIZE-1,"%s",cpl_parameter_get_string(p));
00296     ck0_nomsg(sinfo_qclog_add_string(qclog_tbl,"QC BP-MAP METHOD",key_value,
00297                                       "BP search method","%s"));
00298 
00299 
00300     ck0(sinfo_pro_save_ims(coeffs_imgl,raw,sof,cfg->coeffsCubeName,
00301                PRO_BP_COEFF,qclog_tbl,plugin_id,config),
00302       "cannot dump ims %s", cfg->coeffsCubeName);   
00303 
00304     sinfo_free_table(&qclog_tbl);
00305 
00306     /* =======================================================================
00307        DETERMINES BAD PIXEL MAP
00308        =======================================================================
00309     */
00310     sinfo_msg("Generates bad pixel map");
00311     cknull(mask = sinfo_new_search_bad_pixels (coeffs_imgl, 
00312                             cfg->threshSigmaFactor, 
00313                             cfg->nonlinearThresh, 
00314                             cfg->loReject, 
00315                             cfg->hiReject),
00316        "could not create bad pixel mask!") ;
00317 
00318     sinfo_free_imagelist(&coeffs_imgl) ;
00319 
00320     check_nomsg(n_bad = sinfo_new_count_bad_pixels(mask)) ;
00321     sinfo_msg( "no of bad pixels: %d", n_bad ) ;
00322 
00323     /* QC LOG */
00324     cknull_nomsg(qclog_tbl=sinfo_qclog_init());
00325     check_nomsg(p = cpl_parameterlist_find(config, "sinfoni.bp.method"));
00326     snprintf(key_value, MAX_NAME_SIZE-1,"%s",cpl_parameter_get_string(p));
00327     ck0_nomsg(sinfo_qclog_add_string(qclog_tbl,"QC BP-MAP METHOD",key_value,
00328                                      "BP search method","%s"));
00329 
00330     ck0_nomsg(sinfo_qclog_add_int(qclog_tbl,"QC BP-MAP NBADPIX",n_bad,
00331                                   "No of bad pixels","%d"));
00332 
00333     ck0(sinfo_pro_save_ima(mask,raw,sof,cfg->outName,PRO_BP_MAP_NL,qclog_tbl,
00334                plugin_id,config),
00335     "cannot save ima %s", cfg->outName);
00336 
00337 
00338     /* ======================================================================
00339        FREE MEMORY
00340        ======================================================================
00341     */
00342 
00343 
00344 
00345     sinfo_free_table(&qclog_tbl);
00346     sinfo_free_image(&mask); 
00347     sinfo_free_frameset(&raw); 
00348     sinfo_detlin_free(&cfg);
00349 
00350     return 0;  
00351 
00352  cleanup:
00353     sinfo_free_frameset(&raw_on);
00354     sinfo_free_frameset(&raw_of);
00355 
00356     sinfo_free_image(&mask) ;
00357     sinfoni_free_vector(&cube_mean);
00358     sinfo_free_imagelist(&imgl);
00359     sinfo_free_imagelist(&coeffs_imgl);
00360     sinfo_free_table(&gain);
00361     sinfo_free_table(&qclog_tbl);
00362     sinfo_free_table(&det_lin);
00363     sinfoni_free_vector(&vec_adl);
00364     sinfoni_free_vector(&vec_med);
00365     sinfo_free_polynomial(&pol);
00366     sinfo_detlin_free(&cfg);
00367     sinfo_free_frameset(&raw);
00368 
00369     return -1;
00370 
00371 
00372 }

Generated on Wed Jan 17 08:33:42 2007 for SINFONI Pipeline Reference Manual by  doxygen 1.4.4