bp_lin.c

00001 /*----------------------------------------------------------------------------
00002    
00003    File name    :       bp_lin.c
00004    Author       :   J. Schreiber
00005    Created on   :   May 5, 2003
00006    Description  :   Different methods for searching for bad pixels
00007                         used in the recipe spiffi_badsearch 
00008 
00009  ---------------------------------------------------------------------------*/
00010 
00011 /*----------------------------------------------------------------------------
00012                                 Includes
00013  ---------------------------------------------------------------------------*/
00014 #include "bp_lin.h"
00015 #include "detlin_ini_by_cpl.h"
00016 #include "sinfoni_pro_save.h"
00017 #include "sinfoni_functions.h"
00018 /*----------------------------------------------------------------------------
00019                                 Defines
00020  ---------------------------------------------------------------------------*/
00021 
00022 
00023 /*----------------------------------------------------------------------------
00024                              Function Definitions
00025  ---------------------------------------------------------------------------*/
00026 
00027 /*----------------------------------------------------------------------------
00028    Function     :       badSearchLin()
00029    In           :       ini_file: file name of according .ini file
00030    Out          :       integer (0 if it worked, -1 if it doesn't) 
00031    Job          :       
00032 
00033    this function searches for static bad pixels in stacks of flatfield frames
00034    with in/decreasing intensities. For each pixel position a curve is plotted 
00035    of the pixel intensity in each plane against the clean mean of the whole 
00036    plane.
00037 
00038    A polynomial is fit and the found coefficients are stored sequentially in
00039    a data cube. Then the deviation of the linear coefficients from the clean 
00040    mean is computed for each pixel and the pixel are declared bad if the 
00041    deviations exceed a threshold defined by a factor of the clean standard 
00042    deviation.
00043 
00044    For the resting good pixels the non-linear coefficients are examined.
00045    If one coefficients deviate more than a user given threshold the pixel
00046    is also declared as bad. The data cubus with the fit results and
00047    a bad pixel mask image is stored
00048 
00049  ---------------------------------------------------------------------------*/
00050 int badSearchLin (const char* plugin_id, cpl_parameterlist* config, cpl_frameset* sof)
00051 {
00052     const char* _id = "badSearchLin";
00053     detlin_config * cfg=NULL;
00054 
00055     OneImage ** imagelist ;
00056     OneImage * mask=NULL;
00057     OneImage * eimg=NULL;
00058 
00059     OneCube * cubeflat=NULL;
00060     OneCube * coeffsCube=NULL;
00061 
00062     cpl_image* bp_img=NULL;
00063     cpl_image* bp_imgw=NULL;
00064     cpl_image* cimg=NULL;
00065     cpl_image* wimg=NULL;
00066 
00067     cpl_imagelist* cub_ims=NULL;
00068 
00069     cpl_vector* cube_mean=NULL;
00070     cpl_vector* cube_median=NULL;
00071 
00072     cpl_frameset* raw=NULL;
00073 
00074     cpl_table* qclog_tbl=NULL;
00075 
00076     cpl_parameter *p=NULL;             
00077     cpl_polynomial *pol=NULL;
00078     cpl_vector* vec_adl=NULL;
00079     cpl_vector* vec_med=NULL;
00080     double* mse=NULL;
00081     int i=0;
00082     int n_bad=0 ;
00083     char* key_value=NULL;
00084     char* key_name=NULL;
00085     cpl_table* det_lin=NULL;
00086     int nraw=0;
00087     int* status=NULL;
00088     cpl_table* gain=NULL;
00089     int ngain=0;
00090  
00091     /* parse the file names and parameters to the bad_config data structure cfg */
00092     raw=cpl_frameset_new();
00093     cfg = parse_cpl_input_detlin(config,sof,&raw) ;
00094 
00095 
00096     if (cfg == NULL)
00097     {
00098         cpl_msg_error (_id,"badSearchLin: could not parse .ini file!\n") ;
00099         cpl_frameset_delete(raw);
00100         return -1 ;
00101     }
00102 
00103 
00104     /* ======================================================================
00105        DETERMINES LINEARITY COEFF AS DFO DOES
00106        ======================================================================
00107     */
00108     det_lin=get_linearity(raw);
00109 
00110     nraw=cpl_table_get_nrow(det_lin);
00111     vec_adl=cpl_vector_new(nraw);
00112     vec_med=cpl_vector_new(nraw);
00113 
00114     for(i=0;i<nraw;i++) {
00115       cpl_vector_set(vec_adl,i,cpl_table_get_double(det_lin,"adl",i,status)); 
00116       cpl_vector_set(vec_med,i,cpl_table_get_double(det_lin,"med",i,status));
00117     }
00118     pol=cpl_polynomial_fit_1d_create(vec_adl,vec_med,cfg->order,mse);
00119      cpl_vector_delete(vec_adl);
00120     cpl_vector_delete(vec_med);
00121 
00122     qclog_tbl=sinfoni_qclog_init(cfg->order+2);
00123     key_value = cpl_calloc(FILE_NAME_SZ,sizeof(char));
00124     key_name = cpl_calloc(FILE_NAME_SZ,sizeof(char));
00125 
00126     for(i=0;i<cfg->order+1;i++) {
00127 
00128        sprintf(key_name,"%s%i%s","QC BP-MAP LIN",i," MED");
00129        sprintf(key_value,"%g",cpl_polynomial_get_coeff(pol,&i));
00130        sinfoni_qclog_add(qclog_tbl,i,key_name,"CPL_TYPE_DOUBLE",key_value,
00131                                             "Linearity Polynomial Coeff");
00132     }
00133     cpl_polynomial_delete(pol);
00134     p = cpl_parameterlist_find(config, "sinfoni.bp.method");
00135     sprintf(key_value, "%s",cpl_parameter_get_string(p));
00136     sinfoni_qclog_add(qclog_tbl,cfg->order+1,"QC BP-MAP METHOD","CPL_TYPE_STRING",key_value,
00137                                             "BP search method");
00138 
00139     if(-1 == sinfoni_pro_save_tbl(det_lin,raw,sof,
00140             BP_LIN_LIN_DET_INFO_OUT_FILENAME,"LIN_DET_INFO",qclog_tbl,
00141             plugin_id,config)) {
00142       cpl_msg_error(_id,"cannot dump ims %s",BP_LIN_LIN_DET_INFO_OUT_FILENAME);
00143 
00144       cpl_table_delete(qclog_tbl);
00145       cpl_free(key_name);
00146       cpl_free(key_value);
00147       cpl_table_delete(det_lin);
00148       cpl_frameset_delete(raw);
00149       detlin_free(cfg);
00150       return -1;
00151 
00152     }
00153     cpl_table_delete(det_lin);
00154     cpl_table_delete(qclog_tbl);
00155     cpl_free(key_name);
00156     cpl_free(key_value);
00157 
00158 
00159 
00160     /*=======================================================*/
00161     
00162     cpl_msg_info(_id,"computes gain");
00163    
00164     gain=get_gain3(raw);
00165      
00166     if(gain==NULL) {
00167       cpl_msg_error(_id,"Gain computation is failed");
00168       cpl_frameset_delete(raw);
00169       detlin_free(cfg);
00170       return -1;
00171     }
00172     ngain=cpl_table_get_nrow(gain); 
00173     qclog_tbl=sinfoni_qclog_init(ngain);
00174     key_value = cpl_calloc(FILE_NAME_SZ,sizeof(char));
00175     key_name = cpl_calloc(FILE_NAME_SZ,sizeof(char));
00176     for(i=0;i<ngain;i++) {
00177         sprintf(key_name,"%s%i","QC GAIN",i);
00178         sprintf(key_value,"%g",cpl_table_get_double(gain,"gain",i,status));
00179         sinfoni_qclog_add(qclog_tbl,i,"QC GAIN","CPL_TYPE_DOUBLE",key_value,
00180                                             "Detector gain");
00181     }
00182 
00183    
00184     if(-1 == sinfoni_pro_save_tbl(gain,raw,sof,BP_LIN_GAIN_OUT_FILENAME,
00185             "GAIN_INFO",qclog_tbl,plugin_id,config)) {
00186       cpl_msg_error(_id,"cannot dump tbl %s", BP_LIN_GAIN_OUT_FILENAME);   
00187 
00188       cpl_table_delete(qclog_tbl);
00189       cpl_free(key_name);
00190       cpl_free(key_value);
00191       cpl_table_delete(gain);
00192       cpl_frameset_delete(raw);
00193       detlin_free(cfg);
00194       return -1;
00195 
00196     }
00197     
00198     cpl_table_delete(gain);
00199     cpl_table_delete(qclog_tbl);
00200     cpl_free(key_name);
00201     cpl_free(key_value);
00202    
00203 
00204     /* ===========================================================================
00205        DETERMINES LINEARITY COEFF AS MPE DOES
00206        ===========================================================================
00207     */
00208     imagelist = (OneImage**) cpl_calloc (cfg -> nframes, sizeof(OneImage*)) ;
00209     for ( i = 0 ; i < cfg->nframes ; i++ )
00210     {
00211         if(is_fits_file (cfg->framelist[i]) != 1) {
00212            cpl_msg_error(_id,"Input file %s is not FITS",cfg->framelist[i] );
00213            return -1;
00214         }
00215         imagelist[i] = load_image(cfg->framelist[i]) ;
00216     }
00217 
00218     cubeflat = list_make_cube(imagelist, cfg-> nframes) ;
00219     if ( cubeflat == NULL )
00220     {
00221         cpl_msg_error(_id," could not make a cube out of a list!\n") ;
00222         if(imagelist) cpl_free(imagelist) ;
00223         cpl_frameset_delete(raw);
00224         detlin_free(cfg);
00225     return -1 ;
00226     }
00227  
00228     /*----------------------------------------------------------------
00229      *---------------------- SEARCH FOR BAD PIXELS---------------------
00230      *--------------------------------------------------------------*/
00231     cpl_msg_info(_id,"Search for bad pixels");
00232     coeffsCube = fitIntensityCourse(cubeflat, cfg->order, cfg->loReject, cfg->hiReject) ;
00233    
00234    
00235 
00236     if ( coeffsCube == NULL )
00237     {    
00238     cpl_msg_error (_id,"could not fit polynomial and store coeffs in a data cube!") ;
00239 
00240         destroy_cube (cubeflat) ;
00241         if(imagelist) cpl_free(imagelist) ;
00242         cpl_frameset_delete(raw);
00243         detlin_free(cfg);
00244         return -1 ;
00245     }
00246 
00247     /*---store the polynomial fit coefficients in a data cube----*/
00248     /* CUBE NOT DUMPED BECAUSE PROBLEMS ON LAST PLANE */
00249     cube_mean=cpl_vector_new(coeffsCube->np);
00250     cube_median=cpl_vector_new(coeffsCube->np);
00251     /* cub_ims=cpl_imagelist_new(coeffsCube->lx,coeffsCube->ly,coeffsCube->np,CPL_TYPE_FLOAT); */
00252     cub_ims=cpl_imagelist_new();
00253 
00254 
00255 
00256     /* QC LOG */
00257     qclog_tbl=sinfoni_qclog_init(coeffsCube->np+1);
00258 
00259     key_value = cpl_calloc(FILE_NAME_SZ,sizeof(char));
00260     key_name = cpl_calloc(FILE_NAME_SZ,sizeof(char));
00261 
00262     for(i=0;i<coeffsCube->np;i++) {
00263 
00264        eimg=coeffsCube->plane[i];
00265        wimg=cpl_image_wrap_float(eimg->lx,eimg->ly,eimg->data);
00266        cimg=cpl_image_duplicate(wimg);
00267        cpl_image_unwrap(wimg);
00268        cpl_imagelist_set(cub_ims,cimg,i);
00269        cpl_vector_set(cube_mean,i,cpl_image_get_median(cimg));
00270 
00271        sprintf(key_name,"%s%i%s","QC BP-MAP LIN",i," MEAN");
00272        sprintf(key_value,"%g",cpl_vector_get(cube_mean,i));
00273        sinfoni_qclog_add(qclog_tbl,i,key_name,"CPL_TYPE_DOUBLE",key_value,
00274                                             "Linearity Polynomial Coeff");
00275 
00276     }
00277     p = cpl_parameterlist_find(config, "sinfoni.bp.method");
00278     sprintf(key_value, "%s",cpl_parameter_get_string(p));
00279     sinfoni_qclog_add(qclog_tbl,coeffsCube->np,"QC BP-MAP METHOD","CPL_TYPE_STRING",key_value,
00280                                             "BP search method");
00281 
00282 
00283 
00284     if(-1 == sinfoni_pro_save_ims(cub_ims,raw,sof,cfg->coeffsCubeName,
00285             PRO_BP_COEFF,qclog_tbl,plugin_id,config)) {
00286       cpl_msg_error(_id,"cannot dump ims %s", cfg->coeffsCubeName);   
00287 
00288       cpl_table_delete(qclog_tbl);
00289       cpl_free(key_value);
00290       cpl_free(key_name);
00291        cpl_vector_delete(cube_mean);
00292       cpl_vector_delete(cube_median);
00293       cpl_imagelist_delete(cub_ims); 
00294        cpl_free(imagelist) ;
00295      if(cubeflat)   destroy_cube (cubeflat) ;
00296       if(coeffsCube) destroy_cube(coeffsCube) ;
00297       detlin_free(cfg);
00298       cpl_frameset_delete(raw); 
00299 
00300       return -1;
00301 
00302     }
00303     cpl_table_delete(qclog_tbl);
00304     cpl_free(key_value);
00305     cpl_free(key_name);
00306     
00307 
00308 
00309     /* ===========================================================================
00310        DETERMINES BAD PIXEL MAP
00311        ===========================================================================
00312     */
00313     cpl_msg_info(_id,"Generates bad pixel map");
00314     mask = searchBadPixels (coeffsCube, 
00315                             cfg->threshSigmaFactor, 
00316                             cfg->nonlinearThresh, 
00317                             cfg->loReject, 
00318                             cfg->hiReject) ;
00319     if (mask == NULL)
00320     {    
00321     cpl_msg_error (_id,"could not create bad pixel mask!") ;
00322         cpl_vector_delete(cube_mean);
00323         cpl_vector_delete(cube_median);
00324         cpl_imagelist_delete(cub_ims); 
00325         if(imagelist) cpl_free(imagelist) ;
00326         if(cubeflat)   destroy_cube (cubeflat) ;
00327         if(coeffsCube) destroy_cube(coeffsCube) ;
00328         if(mask)  cpl_free(mask) ;
00329         detlin_free(cfg);
00330         cpl_frameset_delete(raw);
00331         return -1 ;
00332     }
00333 
00334     n_bad = countBadPixels(mask) ;
00335     cpl_msg_info (_id, "no of bad pixels: %d", n_bad ) ;
00336     bp_imgw=cpl_image_wrap_float(mask->lx,mask->ly,mask->data);
00337     bp_img=cpl_image_duplicate(bp_imgw);
00338     cpl_image_unwrap(bp_imgw);
00339 
00340     /* QC LOG */
00341     qclog_tbl=sinfoni_qclog_init(2);
00342     key_value = cpl_calloc(FILE_NAME_SZ,sizeof(char));
00343     p = cpl_parameterlist_find(config, "sinfoni.bp.method");
00344     sprintf(key_value, "%s",cpl_parameter_get_string(p));
00345 
00346     sinfoni_qclog_add(qclog_tbl,0,"QC BP-MAP METHOD","CPL_TYPE_STRING",
00347                        key_value,"BP search method");
00348 
00349     sprintf(key_value,"%d",n_bad);
00350     sinfoni_qclog_add(qclog_tbl,1,"QC BP-MAP NBADPIX","CPL_TYPE_INT",
00351                                   key_value,"No of bad pixels");
00352    
00353 
00354 
00355     if(-1 == sinfoni_pro_save_ima(bp_img,raw,sof,cfg->outName,
00356             PRO_BP_MAP_NL,qclog_tbl,plugin_id,config)) {
00357       cpl_msg_error(_id,"cannot save ima %s", cfg->outName);
00358        cpl_table_delete(qclog_tbl);
00359         cpl_free(key_value);
00360 
00361         if(cubeflat)   destroy_cube(cubeflat) ; 
00362         if(coeffsCube) destroy_cube(coeffsCube) ; 
00363         if(mask)       destroy_image(mask);
00364         cpl_image_delete(bp_img); 
00365     cpl_imagelist_delete(cub_ims);
00366     if(imagelist)  cpl_free(imagelist) ;
00367     cpl_vector_delete(cube_mean);
00368     cpl_vector_delete(cube_median);
00369         cpl_frameset_delete(raw);
00370         detlin_free(cfg);
00371         return -1;
00372 
00373     }
00374 
00375 
00376     /* ===========================================================================
00377        FREE MEMORY
00378        ===========================================================================
00379     */
00380 
00381     cpl_table_delete(qclog_tbl);
00382     cpl_free(key_value);
00383     if(cubeflat)   destroy_cube(cubeflat) ;
00384     if(coeffsCube) destroy_cube(coeffsCube) ;
00385     if(mask)       destroy_image(mask) ;
00386     cpl_imagelist_delete(cub_ims); 
00387     if(imagelist)  cpl_free(imagelist) ;
00388     cpl_image_delete(bp_img); 
00389     cpl_frameset_delete(raw); 
00390     detlin_free(cfg);
00391 
00392     cpl_vector_delete(cube_mean);
00393     cpl_vector_delete(cube_median);
00394 
00395     return 0;  
00396 }
00397 

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