00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 #ifdef HAVE_CONFIG_H
00029 # include <config.h>
00030 #endif
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040 #include "sinfo_utl_cube_combine.h"
00041 #include "sinfo_utilities.h"
00042 #include "sinfo_utilities_scired.h"
00043 #include <stdio.h>
00044 #include "sinfo_key_names.h"
00045 #include "sinfo_error.h"
00046 #include "sinfo_utils_wrappers.h"
00047 #include "sinfo_pro_save.h"
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00070
00077
00078 int sinfo_utl_cube_combine(
00079 cpl_parameterlist * parlist,
00080 cpl_frameset * framelist)
00081 {
00082 cpl_parameter * param=NULL ;
00083 const char * name_o=NULL ;
00084 const char * name_i=NULL ;
00085 const char * name_m=NULL ;
00086
00087 cpl_propertylist * plist=NULL ;
00088 cpl_image * image=NULL ;
00089 cpl_frame * product_frame=NULL;
00090 cpl_frame * prod_frm=NULL;
00091 cpl_frame * frame=NULL;
00092
00093 int z=0;
00094 int z_min=0;
00095 int z_max=0;
00096 int z_siz=0;
00097 int z_stp=100;
00098 int xsize=0;
00099 int ysize=0;
00100 int scales_sky=0;
00101
00102 int i=0;
00103 int j=0;
00104 int n=0;
00105
00106 int size_x=0;
00107 int size_y=0;
00108 cpl_image* j_img=NULL;
00109 cpl_image* m_img=NULL;
00110
00111 cpl_imagelist ** cube_object=NULL;
00112 cpl_imagelist * cube_jitter=NULL;
00113 cpl_imagelist * cube_mask=NULL;
00114
00115 char** files=NULL;
00116 int nframes=0;
00117
00118 double * times=NULL;
00119 float * offsetx=NULL;
00120 float * offsety=NULL;
00121 float ref_offx=0;
00122 float ref_offy=0;
00123 float tmpoffx=0;
00124 float tmpoffy=0;
00125 double kappa=2;
00126 int ks_clip=0;
00127
00128
00129
00130 FILE* file_list=NULL;
00131 int cnt=0;
00132 int onp=0;
00133
00134
00135
00136 check_nomsg(param = cpl_parameterlist_find(parlist,
00137 "sinfoni.sinfo_utl_cube_combine.ks_clip"));
00138 check_nomsg(ks_clip = cpl_parameter_get_bool(param));
00139
00140 check_nomsg(param = cpl_parameterlist_find(parlist,
00141 "sinfoni.sinfo_utl_cube_combine.scale_sky"));
00142 check_nomsg(scales_sky = cpl_parameter_get_bool(param));
00143
00144
00145 check_nomsg(param = cpl_parameterlist_find(parlist,
00146 "sinfoni.sinfo_utl_cube_combine.kappa"));
00147 check_nomsg(kappa = cpl_parameter_get_double(param));
00148
00149 check_nomsg(param = cpl_parameterlist_find(parlist,
00150 "sinfoni.sinfo_utl_cube_combine.name_i"));
00151 check_nomsg(name_i = cpl_parameter_get_string(param));
00152
00153 check_nomsg(param = cpl_parameterlist_find(parlist,
00154 "sinfoni.sinfo_utl_cube_combine.name_o"));
00155 check_nomsg(name_o = cpl_parameter_get_string(param));
00156
00157
00158
00159 check_nomsg(param = cpl_parameterlist_find(parlist,
00160 "sinfoni.sinfo_utl_cube_combine.xsize"));
00161 check_nomsg(xsize = cpl_parameter_get_int(param));
00162
00163 check_nomsg(param = cpl_parameterlist_find(parlist,
00164 "sinfoni.sinfo_utl_cube_combine.ysize"));
00165 check_nomsg(ysize = cpl_parameter_get_int(param));
00166
00167
00168
00169
00170 if ( NULL == (file_list = fopen (name_i, "r" ) ) )
00171 {
00172 sinfo_msg_error("cannot open %s\n", name_i) ;
00173 goto cleanup ;
00174 }
00175
00176 cnt = 0 ;
00177 while ( fscanf( file_list, "%f %f",&tmpoffx, &tmpoffy) != EOF )
00178 {
00179 cnt ++ ;
00180 }
00181 fclose(file_list);
00182
00183 nframes= cnt ;
00184 cknull(times = (double*) cpl_calloc (nframes, sizeof (double)),
00185 " could not allocate memory!") ;
00186
00187 cknull(offsetx = (float*) cpl_calloc (nframes, sizeof(float)),
00188 " could not allocate memory!") ;
00189
00190 cknull(offsety = (float*) cpl_calloc (nframes, sizeof(float)),
00191 " could not allocate memory!") ;
00192
00193 files = (char**) cpl_calloc(MAX_NAME_SIZE,sizeof(char*));
00194 if ( NULL == (file_list = fopen (name_i, "r" ) ) )
00195 {
00196 sinfo_msg_error("cannot open %s\n", name_i) ;
00197 return -1 ;
00198 }
00199
00200 cnt=0;
00201 while ( fscanf( file_list, "%f %f",&tmpoffx,&tmpoffy ) != EOF ) {
00202 frame=cpl_frameset_get_frame(framelist,cnt);
00203 files[cnt]=(char*) cpl_frame_get_filename(frame);
00204 offsetx[cnt]=tmpoffx;
00205 offsety[cnt]=tmpoffy;
00206 cnt ++ ;
00207 }
00208 nframes=cnt;
00209 fclose(file_list);
00210
00211 ck0(sinfo_auto_size_cube3(files,nframes,&ref_offx,&ref_offy,&size_x,&size_y),
00212 "Error resizing cube");
00213 sinfo_msg("output ima size=%dx%d",size_x,size_y);
00214 cknull(cube_object = cpl_calloc(nframes,sizeof(cpl_imagelist*)),
00215 "could not allocate memory") ;
00216
00217 for(n=0;n<nframes;n++) {
00218 cube_object[i]=cpl_imagelist_new();
00219 }
00220
00221 for (i=0;i<nframes;i++) {
00222 times[i]=sinfo_pfits_get_exptime(files[i]);
00223 cube_object[i]=cpl_imagelist_load(files[i],CPL_TYPE_FLOAT,0);
00224 ck0(sinfo_new_assign_offset(i,files[i],offsetx,
00225 offsety,ref_offx,ref_offy),
00226 "Error assigning offsets");
00227 }
00228 onp = cpl_imagelist_get_size(cube_object[0]);
00229
00230 check(cube_jitter = cpl_imagelist_new(),"allocating new data cube object");
00231 check(cube_mask = cpl_imagelist_new(),"allocating new data cube mask");
00232
00233 check(plist=cpl_propertylist_load(files[0],0),
00234 "Cannot read the FITS header") ;
00235
00236 if(scales_sky == 1) {
00237 sinfo_msg("Subtract spatial sinfo_median to each cube plane");
00238 for(n=0;n<nframes;n++) {
00239 sinfo_msg("process cube %d\n",n);
00240 sinfo_new_sinfoni_correct_median_it(&(cube_object[n]));
00241 }
00242 }
00243
00244 for(z=0;z<onp;z+=z_stp) {
00245 z_siz=(z_stp < (onp-z)) ? z_stp : (onp-z);
00246 z_min=z;
00247 z_max=z_min+z_siz;
00248 sinfo_msg("Coadding cubes: range [%4.4d,%4.4d) of 0-%d\n",
00249 z_min,z_max,onp);
00250
00251 for(j=z_min;j<z_max;j++) {
00252
00253 check_nomsg(j_img=cpl_image_new(size_x,size_y,CPL_TYPE_FLOAT));
00254 check_nomsg(cpl_imagelist_set(cube_jitter,j_img,j));
00255 check_nomsg(m_img = cpl_image_new(size_x,size_y,CPL_TYPE_FLOAT));
00256 check_nomsg(cpl_imagelist_set(cube_mask,m_img,j));
00257 }
00258 if(ks_clip == 1){
00259 sinfo_new_combine_jittered_cubes_thomas_range(cube_object,
00260 cube_jitter,
00261 cube_mask,
00262 nframes,
00263 offsetx,offsety,
00264 times,
00265 (char*) "tanh",
00266 z_min,
00267 z_max,
00268 kappa);
00269 } else {
00270 sinfo_new_combine_jittered_cubes_range(cube_object,
00271 cube_jitter,
00272 cube_mask,
00273 nframes,
00274 offsetx,
00275 offsety,
00276 times,
00277 (char*) "tanh",
00278 z_min,
00279 z_max) ;
00280 }
00281
00282 }
00283 sinfo_new_convert_0_to_ZERO_for_cubes(cube_jitter) ;
00284
00285 name_m="out_cube_mask.fits";
00286
00287
00288
00289 check_nomsg(product_frame = cpl_frame_new());
00290 check_nomsg(cpl_frame_set_filename(product_frame, name_o)) ;
00291 check_nomsg(cpl_frame_set_tag(product_frame, SI_UTL_CUBE_COMBINE_PROCUBE));
00292 check_nomsg(cpl_frame_set_type(product_frame, CPL_FRAME_TYPE_IMAGE)) ;
00293 check_nomsg(cpl_frame_set_group(product_frame, CPL_FRAME_GROUP_PRODUCT)) ;
00294 check_nomsg(cpl_frame_set_level(product_frame, CPL_FRAME_LEVEL_FINAL)) ;
00295
00296
00297
00298
00299
00300
00301
00302
00303
00304 check(cpl_imagelist_save(cube_jitter, name_o, CPL_BPP_DEFAULT, plist,
00305 CPL_IO_DEFAULT),"Could not save product");
00306 check_nomsg(cpl_frameset_insert(framelist, product_frame));
00307
00308
00309 check_nomsg(prod_frm = cpl_frame_new());
00310 check_nomsg(cpl_frame_set_filename(prod_frm, name_m)) ;
00311 check_nomsg(cpl_frame_set_tag(prod_frm, SI_UTL_CUBE_COMBINE_PROMASK)) ;
00312 check_nomsg(cpl_frame_set_type(prod_frm, CPL_FRAME_TYPE_IMAGE)) ;
00313 check_nomsg(cpl_frame_set_group(prod_frm, CPL_FRAME_GROUP_PRODUCT)) ;
00314 check_nomsg(cpl_frame_set_level(prod_frm, CPL_FRAME_LEVEL_FINAL)) ;
00315
00316 check(cpl_imagelist_save(cube_mask, name_m, CPL_BPP_DEFAULT, plist,
00317 CPL_IO_DEFAULT),"Could not save product");
00318 check_nomsg(cpl_frameset_insert(framelist, prod_frm));
00319
00320 cleanup:
00321 sinfo_free_imagelist(&cube_jitter) ;
00322 sinfo_free_image(&image) ;
00323 sinfo_free_imagelist(&cube_mask) ;
00324
00325 if(cube_object != NULL) {
00326 for (i=0;i< nframes;i++){
00327 sinfo_free_imagelist(&(cube_object[i]));
00328 }
00329 sinfo_free_array_imagelist(&cube_object);
00330 }
00331 sinfo_free_propertylist(&plist);
00332 if(files != NULL) {
00333
00334
00335
00336
00337
00338
00339
00340
00341 cpl_free(files);
00342 files=NULL;
00343 }
00344 sinfo_free_double(×);
00345 sinfo_free_float(&offsetx);
00346 sinfo_free_float(&offsety);
00347
00348
00349 if (cpl_error_get_code()) {
00350 return -1 ;
00351 }
00352 else {
00353 return 0 ;
00354 }
00355 }