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
00029
00030
00031
00032
00033
00034 #include "utl_cube_combine.h"
00035 #include "eclipse.h"
00036 #include "utilities.h"
00037 #include <stdio.h>
00038 #include "sinfoni_key_names.h"
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00060
00061 int si_utl_cube_combine(
00062 cpl_parameterlist * parlist,
00063 cpl_frameset * framelist)
00064 {
00065 const char * fctid = "si_utl_cube_combine" ;
00066 cpl_parameter * param=NULL ;
00067 const char * operation=NULL ;
00068 const char * method=NULL ;
00069 double sigma=0 ;
00070 cpl_image * ima1=NULL ;
00071 cpl_image * ima2=NULL ;
00072
00073 const char * name_o=NULL ;
00074 const char * name_i=NULL ;
00075 const char * name_m=NULL ;
00076
00077 cpl_propertylist * plist=NULL ;
00078 cpl_image * image=NULL ;
00079 cpl_frame * product_frame=NULL;
00080 cpl_frame * frame=NULL;
00081
00082
00083 int xsize=0;
00084 int ysize=0;
00085
00086 int i=0;
00087 int n=0;
00088
00089 OneCube ** cube_object;
00090 OneCube * cube_jitter=NULL;
00091 OneCube * cube_mask=NULL;
00092 OneCube ** cube_out;
00093 char** files;
00094 int nframes=0;
00095 float time=0;
00096 float * times=NULL;
00097 float * offsetx=NULL;
00098 float * offsety=NULL;
00099 float * offs_x=NULL;
00100 float * offs_y=NULL;
00101
00102 char ** fileoffs;
00103 char * filename=NULL;
00104 char * tmpname=NULL;
00105 char * tmptag=NULL;
00106 float tmpoffx=0;
00107 float tmpoffy=0;
00108
00109 cpl_imagelist * cub_jit=NULL;
00110 cpl_imagelist * cub_msk=NULL;
00111 OneImage* eimg=NULL;
00112 cpl_image * cimg=NULL;
00113
00114 FILE* file_list=NULL;
00115 int cnt=0;
00116 int np=0;
00117
00118
00119
00120 param = cpl_parameterlist_find(parlist, "sinfoni.si_utl_cube_combine.op");
00121 operation = cpl_parameter_get_string(param);
00122
00123 param = cpl_parameterlist_find(parlist, "sinfoni.si_utl_cube_combine.method");
00124 method = cpl_parameter_get_string(param);
00125
00126
00127 param = cpl_parameterlist_find(parlist, "sinfoni.si_utl_cube_combine.name_o");
00128 name_o = cpl_parameter_get_string(param);
00129
00130
00131 param = cpl_parameterlist_find(parlist, "sinfoni.si_utl_cube_combine.name_i");
00132 name_i = cpl_parameter_get_string(param);
00133
00134
00135
00136 param = cpl_parameterlist_find(parlist,"sinfoni.si_utl_cube_combine.xsize");
00137 xsize = cpl_parameter_get_int(param) ;
00138
00139
00140 param = cpl_parameterlist_find(parlist,"sinfoni.si_utl_cube_combine.ysize");
00141 ysize = cpl_parameter_get_int(param) ;
00142
00143
00144 param = cpl_parameterlist_find(parlist,"sinfoni.si_utl_cube_combine.sigma");
00145 sigma = cpl_parameter_get_double(param) ;
00146
00147
00148
00149
00150
00151 if ( NULL == (file_list = fopen (name_i, "r" ) ) )
00152 {
00153 cpl_msg_error(fctid,"cannot open %s\n", name_i) ;
00154 return -1 ;
00155 }
00156 cnt = 0 ;
00157 tmpname = (char*) cpl_calloc(MAX_NAME_SIZE,sizeof(char));
00158 tmptag = (char*) cpl_calloc(MAX_NAME_SIZE,sizeof(char));
00159 while ( fscanf( file_list, "%f %f",&tmpoffx, &tmpoffy) != EOF )
00160 {
00161 cnt ++ ;
00162 }
00163
00164
00165 fclose(file_list);
00166 n= cnt ;
00167 files = (char**) cpl_calloc(MAX_NAME_SIZE,sizeof(char*));
00168 offs_x = (float*) cpl_calloc(cnt,sizeof(float*));
00169 offs_y = (float*) cpl_calloc(cnt,sizeof(float*));
00170 cnt=0;
00171 if ( NULL == (file_list = fopen (name_i, "r" ) ) )
00172 {
00173 cpl_msg_error(fctid,"cannot open %s\n", name_i) ;
00174 return -1 ;
00175 }
00176
00177 while ( fscanf( file_list, "%f %f",&tmpoffx,&tmpoffy ) != EOF ) {
00178 frame=cpl_frameset_get_frame(framelist,cnt);
00179 files[cnt]=cpl_frame_get_filename(frame);
00180 offs_x[cnt]=tmpoffx;
00181 offs_x[cnt]=tmpoffy;
00182 printf("filename=%s\n",files[cnt]);
00183 cnt ++ ;
00184 }
00185 nframes=cnt;
00186 fclose(file_list);
00187
00188
00189 times = new_float_array(nframes);
00190 offsetx = new_float_array(nframes);
00191 offsety = new_float_array(nframes);
00192
00193 fileoffs = (char**) cpl_calloc(nframes,sizeof(char*));
00194 filename = (char*) cpl_calloc(MAX_NAME_SIZE,sizeof(char));
00195 for (n=0;n<nframes;n++){
00196 array_set_value(offsetx, offs_x[n], n);
00197 array_set_value(offsety, offs_y[n], n);
00198 }
00199
00200
00201
00202
00203 printf("ok1\n");
00204
00205 cube_object = new_cube_list(nframes);
00206 if(cube_object == NULL) {
00207 cpl_msg_error(fctid,"could not allocate memory") ;
00208 cpl_free(fileoffs);
00209 cpl_free(filename);
00210 destroy_array(times);
00211 destroy_array(offsetx);
00212 destroy_array(offsety);
00213 cpl_free(files);
00214 cpl_free(offs_x);
00215 cpl_free(offs_y);
00216 cpl_free(tmpname);
00217 cpl_free(tmptag);
00218 return -1 ;
00219 }
00220
00221 cube_out = new_cube_list(nframes);
00222 for (i=0;i<nframes;i++) {
00223 printf("file=%s\n",files[i]);
00224 time = spiffi_get_exptime(files[i]);
00225 array_set_value(times, time, i);
00226 cube_out[i]=load_cube(files[i]);
00227 insert_cube(cube_object, cube_out[i], i);
00228 }
00229 np = cube_out[0]->np;
00230
00231
00232 for (i=0;i< nframes;i++){
00233 destroy_cube(cube_out[i]);
00234 }
00235 destroy_cube_list(cube_out);
00236 printf("ok2\n");
00237
00238
00239 cube_jitter = newCube(xsize, ysize, np);
00240 if (cube_jitter == NULL){
00241 cpl_msg_error (fctid,"could allocate new data cube");
00242
00243
00244 destroy_cube_list(cube_object);
00245 destroy_cube(cube_jitter);
00246 cpl_free(fileoffs);
00247 cpl_free(filename);
00248 destroy_array(times);
00249 destroy_array(offsetx);
00250 destroy_array(offsety);
00251 cpl_free(files);
00252 cpl_free(offs_x);
00253 cpl_free(offs_y);
00254 cpl_free(tmpname);
00255 cpl_free(tmptag);
00256 return -1;
00257 }
00258 printf("ok3\n");
00259 printf("file0=%s\n",files[0]);
00260 if ((plist=cpl_propertylist_load(files[0],0)) == NULL) {
00261 cpl_msg_error(fctid, "Cannot read the FITS header") ;
00262 return -1 ;
00263 }
00264 if(strcmp(method,"clean_mean")==0) {
00265 cube_mask = combineJitteredCubes(cube_object, cube_jitter, nframes, offsetx,
00266 offsety, times, (char*) "tanh");
00267 } else {
00268 cube_mask = combineCubes(cube_object, cube_jitter, nframes,
00269 offsetx, offsety, sigma, (char*) "tanh");
00270 }
00271 printf("ok4\n");
00272 if (cube_mask == NULL){
00273 cpl_msg_error (fctid,"could not merge the jittered data cubes\n");
00274 return -1;
00275 }
00276 name_m="out_cube_mask.fits";
00277 printf("ok5\n");
00278 save_cube_to_fits_hdr_dump( cube_jitter, (char*)name_o, NULL );
00279 save_cube_to_fits_hdr_dump( cube_mask, (char*) name_m, NULL);
00280 printf("ok6\n");
00281
00282
00283
00284
00285
00286 product_frame = cpl_frame_new();
00287 cpl_frame_set_filename(product_frame, name_o) ;
00288 cpl_frame_set_tag(product_frame, SI_UTL_CUBE_COMBINE_PROCUBE) ;
00289 cpl_frame_set_type(product_frame, CPL_FRAME_TYPE_IMAGE) ;
00290 cpl_frame_set_group(product_frame, CPL_FRAME_GROUP_PRODUCT) ;
00291 cpl_frame_set_level(product_frame, CPL_FRAME_LEVEL_FINAL) ;
00292 if (cpl_error_get_code()) {
00293 cpl_msg_error(fctid, "Error while initialising the product frame") ;
00294 cpl_propertylist_delete(plist) ;
00295 cpl_frame_delete(product_frame) ;
00296 cpl_image_delete(image) ;
00297 return -1 ;
00298 }
00299
00300 printf("ok7\n");
00301
00302
00303
00304
00305
00306
00307
00308
00309
00310
00311
00312
00313
00314
00315
00316
00317
00318 cub_jit=cpl_imagelist_new(cube_jitter->lx,cube_jitter->ly,cube_jitter->np,CPL_TYPE_FLOAT);
00319 for(i=0;i<cube_jitter->np;i++) {
00320 eimg=cube_jitter->plane[i];
00321 cimg=cpl_image_wrap_float(eimg->lx,eimg->ly,eimg->data);
00322 cpl_imagelist_set(cub_jit,cimg,i);
00323 }
00324
00325
00326 if (cpl_imagelist_save(cub_jit, name_o, CPL_BPP_DEFAULT, plist,
00327 CPL_IO_DEFAULT) != CPL_ERROR_NONE) {
00328 cpl_msg_error(fctid, "Could not save product");
00329 cpl_propertylist_delete(plist) ;
00330 cpl_frame_delete(product_frame) ;
00331 cpl_imagelist_delete(cub_jit) ;
00332 return -1 ;
00333 }
00334 cpl_frameset_insert(framelist, product_frame) ;
00335
00336
00337
00338
00339 product_frame = cpl_frame_new();
00340 cpl_frame_set_filename(product_frame, name_m) ;
00341 cpl_frame_set_tag(product_frame, SI_UTL_CUBE_COMBINE_PROMASK) ;
00342 cpl_frame_set_type(product_frame, CPL_FRAME_TYPE_IMAGE) ;
00343 cpl_frame_set_group(product_frame, CPL_FRAME_GROUP_PRODUCT) ;
00344 cpl_frame_set_level(product_frame, CPL_FRAME_LEVEL_FINAL) ;
00345 if (cpl_error_get_code()) {
00346 cpl_msg_error(fctid, "Error while initialising the product frame") ;
00347 cpl_propertylist_delete(plist) ;
00348 cpl_frame_delete(product_frame) ;
00349 cpl_image_delete(image) ;
00350 return -1 ;
00351 }
00352
00353 cub_msk=cpl_imagelist_new(cube_mask->lx,cube_mask->ly,cube_mask->np,CPL_TYPE_FLOAT);
00354 for(i=0;i<cube_mask->np;i++) {
00355 eimg=cube_mask->plane[i];
00356 cimg=cpl_image_wrap_float(eimg->lx,eimg->ly,eimg->data);
00357 cpl_imagelist_set(cub_msk,cimg,i);
00358 }
00359
00360
00361 if (cpl_imagelist_save(cub_msk, name_m, CPL_BPP_DEFAULT, plist,
00362 CPL_IO_DEFAULT) != CPL_ERROR_NONE) {
00363 cpl_msg_error(fctid, "Could not save product");
00364 cpl_propertylist_delete(plist) ;
00365 cpl_frame_delete(product_frame) ;
00366 cpl_imagelist_delete(cub_msk) ;
00367 return -1 ;
00368 }
00369 cpl_frameset_insert(framelist, product_frame) ;
00370
00371
00372
00373 cpl_propertylist_delete(plist) ;
00374 cpl_image_delete(ima1);
00375 cpl_image_delete(ima2);
00376 cpl_image_delete(image) ;
00377
00378
00379
00380 destroy_array(times);
00381 destroy_array(offsetx);
00382 destroy_array(offsety);
00383 destroy_cube(cube_mask);
00384 destroy_cube(cube_jitter);
00385 destroy_cube_list(cube_object);
00386
00387
00388
00389
00390 cpl_free(fileoffs);
00391 cpl_free(filename);
00392 cpl_free(offs_x);
00393 cpl_free(offs_y);
00394
00395
00396
00397 if (cpl_error_get_code()) {
00398 return -1 ;
00399 }
00400 else {
00401 return 0 ;
00402 }
00403 }