00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019 #ifdef HAVE_CONFIG_H
00020 # include <config.h>
00021 #endif
00022 #include "sinfo_utilities.h"
00023 #include "sinfo_error.h"
00024 #include "sinfo_utils_wrappers.h"
00025 #include "sinfo_function_1d.h"
00037 int
00038 sinfo_table_column_dump(cpl_table* t, const char* name, cpl_type type)
00039 {
00040 int nrow=0;
00041 int i=0;
00042 int* pi=NULL;
00043 float* pf=NULL;
00044 double* pd=NULL;
00045 char** ps=NULL;
00046 int index=0;
00047
00048 nrow=cpl_table_get_nrow(t);
00049
00050 switch(type) {
00051
00052 case CPL_TYPE_INT:
00053 pi=cpl_table_get_data_int(t,name);
00054 for(i=0;i<nrow;i++) {
00055 sinfo_msg("val=%d",pi[i]);
00056 }
00057 break;
00058 case CPL_TYPE_FLOAT:
00059 pf=cpl_table_get_data_float(t,name);
00060 for(i=0;i<nrow;i++) {
00061 sinfo_msg("val=%g",pf[i]);
00062 }
00063 break;
00064 case CPL_TYPE_DOUBLE:
00065 pd=cpl_table_get_data_double(t,name);
00066 for(i=0;i<nrow;i++) {
00067 sinfo_msg("val=%g",pd[i]);
00068 }
00069 break;
00070 case CPL_TYPE_STRING:
00071 ps=cpl_table_get_data_string(t,name);
00072 for(i=0;i<nrow;i++) {
00073 sinfo_msg("val=%s",ps[i]);
00074 }
00075 break;
00076 default:
00077 sinfo_msg_error("Wrong column type");
00078 cpl_error_set(cpl_func, CPL_ERROR_TYPE_MISMATCH);
00079 return index;
00080
00081 }
00082 return 0;
00083 }
00084
00085
00086
00096 cpl_table*
00097 sinfo_table_shift_column_spline3(cpl_table* t,
00098 const char* col,
00099 const double shift)
00100 {
00101 cpl_table* out=NULL;
00102 int nrow=0;
00103 int i=0;
00104 int z=0;
00105
00106 float sum=0;
00107 float new_sum=0;
00108
00109 float* pi=NULL;
00110 float* po=NULL;
00111 float* eval=NULL;
00112 float* xnum=NULL;
00113 float* spec=NULL;
00114 float* corrected_spec=NULL;
00115
00116 cknull(t,"null input table");
00117 out=cpl_table_duplicate(t);
00118
00119 nrow=cpl_table_get_nrow(t);
00120 check_nomsg(cpl_table_cast_column(t,col,"FINT",CPL_TYPE_FLOAT));
00121 check_nomsg(cpl_table_cast_column(out,col,"FINT",CPL_TYPE_FLOAT));
00122 pi=cpl_table_get_data_float(t,"FINT");
00123 po=cpl_table_get_data_float(out,"FINT");
00124
00125
00126
00127 xnum=cpl_calloc(nrow,sizeof(float)) ;
00128
00129 for ( i = 0 ; i < nrow ; i++ ) {
00130 xnum[i] = i ;
00131 }
00132
00133 spec=cpl_calloc(nrow,sizeof(float)) ;
00134 corrected_spec=cpl_calloc(nrow,sizeof(float)) ;
00135 eval=cpl_calloc(nrow,sizeof(float)) ;
00136
00137 sum = 0. ;
00138 for ( z = 0 ; z < nrow ; z++ ) {
00139 spec[z] = pi[z] ;
00140 if (isnan(spec[z]) ) {
00141 for ( i = z-1 ; i <= z+1 ; i++ ) {
00142 if ( i < 0 ) continue ;
00143 if ( i >= nrow) continue ;
00144 corrected_spec[i] = ZERO ;
00145 }
00146 spec[z] = 0. ;
00147 }
00148 sum += spec[z] ;
00149 eval[z] = (float)shift+(float)z ;
00150 }
00151
00152 if ( -1 == sinfo_function1d_natural_spline(xnum,spec, nrow,
00153 eval,corrected_spec, nrow))
00154 {
00155 sinfo_msg_error("error in spline interpolation!") ;
00156 goto cleanup;
00157 }
00158
00159 new_sum = 0. ;
00160 for ( z = 0 ; z < nrow ; z++ ) {
00161 if ( isnan(corrected_spec[z]) ) {
00162 continue ;
00163 }
00164 new_sum += corrected_spec[z] ;
00165 }
00166
00167 for ( z = 0 ; z < nrow ; z++ ) {
00168 if ( new_sum == 0. ) new_sum =1. ;
00169 {
00170 if ( isnan(corrected_spec[z]) ) {
00171 po[z] = ZERO ;
00172 } else {
00173 corrected_spec[z] *= sum / new_sum ;
00174 po[z] = corrected_spec[z] ;
00175 }
00176 }
00177 }
00178
00179 sinfo_free_float(&xnum);
00180 sinfo_free_float(&spec) ;
00181 sinfo_free_float(&corrected_spec) ;
00182 sinfo_free_float(&eval) ;
00183
00184 check_nomsg(cpl_table_erase_column(t,"FINT"));
00185 check_nomsg(cpl_table_erase_column(out,col));
00186 check_nomsg(cpl_table_cast_column(out,"FINT",col,CPL_TYPE_DOUBLE));
00187 check_nomsg(cpl_table_erase_column(out,"FINT"));
00188
00189 return out;
00190 cleanup:
00191
00192 sinfo_free_float(&xnum);
00193 sinfo_free_float(&spec) ;
00194 sinfo_free_float(&corrected_spec) ;
00195 sinfo_free_float(&eval) ;
00196 sinfo_free_table(&out);
00197 return NULL;
00198
00199
00200 }
00201
00202
00212 cpl_table*
00213 sinfo_table_shift_column_int(const cpl_table* t,
00214 const char* col,
00215 const double s,
00216 double* r)
00217 {
00218 cpl_table* out=NULL;
00219 int is=(int)s;
00220 int nrow=0;
00221 int i=0;
00222
00223 double* pi=NULL;
00224 double* po=NULL;
00225
00226 cknull(t,"null input table");
00227 out=cpl_table_duplicate(t);
00228 *r=s-is;
00229 nrow=cpl_table_get_nrow(t);
00230 pi=cpl_table_get_data_double(t,col);
00231 po=cpl_table_get_data_double(out,col);
00232 for(i=0;i<nrow;i++) {
00233 if( ((i-is) >=0) && ((i-is) < nrow)) {
00234 po[i-is]=pi[i];
00235 }
00236 }
00237 return out;
00238 cleanup:
00239 sinfo_free_table(&out);
00240 return NULL;
00241
00242 }
00243
00244
00255 cpl_table*
00256 sinfo_table_shift_column_poly(cpl_table* t,
00257 const char* col,
00258 const double shift,
00259 const int order)
00260 {
00261 cpl_table* out=NULL;
00262 int nrow=0;
00263 int i=0;
00264 int flag=0;
00265 int n_points=0;
00266 int firstpos=0;
00267 int z=0;
00268 float eval=0;
00269 float sum=0;
00270 float new_sum=0;
00271 float* pi=NULL;
00272 float* po=NULL;
00273 float* spec=NULL ;
00274 float* corrected_spec=NULL ;
00275 float* xnum=NULL ;
00276 float* tableptr=NULL;
00277
00278 cknull(t,"null input table");
00279 if ( order <= 0 ) {
00280 sinfo_msg_error("wrong order of interpolation polynom given!") ;
00281 goto cleanup;
00282 }
00283
00284 out=cpl_table_duplicate(t);
00285
00286 nrow=cpl_table_get_nrow(t);
00287 cpl_table_cast_column(t,col,"FINT",CPL_TYPE_FLOAT);
00288 cpl_table_cast_column(out,col,"FINT",CPL_TYPE_FLOAT);
00289 pi=cpl_table_get_data_float(t,"FINT");
00290 po=cpl_table_get_data_float(out,"FINT");
00291
00292 n_points = order + 1 ;
00293 if ( n_points % 2 == 0 ) {
00294 firstpos = (int)(n_points/2) - 1 ;
00295 } else {
00296 firstpos = (int)(n_points/2) ;
00297 }
00298 spec=cpl_calloc(nrow,sizeof(float)) ;
00299 corrected_spec=cpl_calloc(nrow,sizeof(float)) ;
00300 xnum=cpl_calloc(order+1,sizeof(float)) ;
00301
00302 for ( i = 0 ; i < n_points ; i++ ) {
00303 xnum[i] = i ;
00304 }
00305
00306
00307 for(i=0;i<nrow;i++) {
00308 corrected_spec[i] = 0. ;
00309 }
00310
00311 sum = 0. ;
00312 for ( z = 0 ; z < nrow ; z++ ) {
00313 spec[z] = pi[z] ;
00314 if (isnan(spec[z]) ) {
00315 spec[z] = 0. ;
00316
00317 for ( i = z - firstpos ; i < z-firstpos+n_points ; i++ ) {
00318 if ( i < 0 ) continue ;
00319 if ( i >= nrow) continue ;
00320 corrected_spec[i] = ZERO ;
00321 }
00322 }
00323 if ( z != 0 && z != nrow - 1 ) {
00324 sum += spec[z] ;
00325 }
00326 }
00327
00328 new_sum = 0. ;
00329 for ( z = 0 ; z < nrow ; z++ ) {
00330
00331
00332
00333
00334
00335
00336 if (isnan(corrected_spec[z])) continue ;
00337 if ( z - firstpos < 0 ) {
00338 tableptr = &spec[0] ;
00339 eval = shift + z ;
00340 } else if ( z - firstpos + n_points >= nrow ) {
00341 tableptr = &spec[nrow - n_points] ;
00342 eval = shift + z + n_points - nrow ;
00343 } else {
00344 tableptr = &spec[z-firstpos] ;
00345 eval = shift + firstpos ;
00346 }
00347
00348 flag=0;
00349 corrected_spec[z]=sinfo_new_nev_ille(xnum,tableptr,order,eval,&flag);
00350 if ( z != 0 && z != nrow - 1 ) {
00351 new_sum += corrected_spec[z] ;
00352 }
00353 }
00354
00355
00356 for (z = 0 ; z < nrow ; z++ ) {
00357 if ( new_sum == 0. ) {
00358 new_sum = 1. ;
00359 }
00360 if ( z == 0 ) {
00361 po[z] = ZERO ;
00362 } else if ( z == nrow - 1 ) {
00363 po[z] = ZERO ;
00364 } else if ( isnan(corrected_spec[z]) ) {
00365 po[z] = ZERO ;
00366 } else {
00367 corrected_spec[z] *= sum / new_sum ;
00368 po[z] = corrected_spec[z] ;
00369 }
00370 }
00371 check_nomsg(cpl_table_erase_column(t,"FINT"));
00372 check_nomsg(cpl_table_erase_column(out,col));
00373 check_nomsg(cpl_table_cast_column(out,"FINT",col,CPL_TYPE_DOUBLE));
00374 check_nomsg(cpl_table_erase_column(out,"FINT"));
00375
00376 sinfo_free_float(&spec) ;
00377 sinfo_free_float(&corrected_spec) ;
00378 sinfo_free_float(&xnum) ;
00379
00380 return out;
00381 cleanup:
00382
00383
00384 sinfo_free_float(&spec) ;
00385 sinfo_free_float(&corrected_spec) ;
00386 sinfo_free_float(&xnum) ;
00387 sinfo_free_table(&out);
00388 return NULL;
00389
00390
00391 }
00392
00393
00394
00395
00396 void sinfo_new_array_set_value( float * array, float value, int i )
00397 {
00398 array[i] = value ;
00399 }
00400 float sinfo_new_array_get_value( float * array, int i )
00401 {
00402 return array[i] ;
00403 }
00404
00405
00406
00407 void sinfo_new_destroy_array(float ** array)
00408 {
00409 if(*array != NULL) {
00410 cpl_free( *array ) ;
00411 *array = NULL;
00412 }
00413 }
00414
00415 void sinfo_new_destroy_stringarray(char ** array, int size_x)
00416 {
00417 int i ;
00418
00419 for ( i = 0 ; i < size_x ; i++ )
00420 {
00421 cpl_free( array[i] ) ;
00422 }
00423 cpl_free( array ) ;
00424 }
00425
00426 void sinfo_new_destroy_2Dintarray(int *** array, int size_x)
00427 {
00428 int i ;
00429
00430 if((*array) != NULL) {
00431 for ( i = 0 ; i < size_x ; i++ ) {
00432 if((*array)[i] != NULL) {
00433 cpl_free( (*array)[i] );
00434 (*array)[i]=NULL;
00435 }
00436 }
00437 cpl_free( *array ) ;
00438 *array=NULL;
00439 }
00440
00441 }
00442
00443
00444 void sinfo_new_intarray_set_value( int * array, int value, int i )
00445 {
00446 array[i] = value ;
00447 }
00448 float sinfo_new_array2D_get_value( float ** array, int x, int y )
00449 {
00450 return array[x][y] ;
00451 }
00452 int ** sinfo_new_2Dintarray( int size_x, int size_y)
00453 {
00454 int ** retVal ;
00455 int i ;
00456
00457 retVal = (int **) cpl_calloc( size_x, sizeof (int*) ) ;
00458 for ( i = 0 ; i < size_x ; i++ )
00459 {
00460 retVal[i] = (int *) cpl_calloc( size_y, sizeof (int)) ;
00461 }
00462 return retVal ;
00463 }
00464
00465 float * sinfo_new_floatarray( int size)
00466 {
00467 return (float *) cpl_calloc( size, sizeof (float) ) ;
00468 }
00469
00470
00471 void sinfo_new_destroy_2Dfloatarray(float *** array, int size_x)
00472 {
00473 int i ;
00474 if((*array) != NULL) {
00475 for ( i = 0 ; i < size_x ; i++ ) {
00476 if((*array)[i] != NULL) {
00477 cpl_free( (*array)[i] );
00478 (*array)[i]=NULL;
00479 }
00480 }
00481 cpl_free( *array ) ;
00482 *array=NULL;
00483 }
00484 }
00485
00486 void sinfo_new_destroy_2Ddoublearray(double *** array, int size_x)
00487 {
00488 int i ;
00489
00490 if((*array) != NULL) {
00491 for ( i = 0 ; i < size_x ; i++ ) {
00492 if((*array)[i] != NULL) {
00493 cpl_free( (*array)[i] );
00494 (*array)[i]=NULL;
00495 }
00496 }
00497 cpl_free( *array ) ;
00498 *array=NULL;
00499 }
00500
00501 }
00502
00503
00504 void sinfo_new_array2D_set_value( float ** array, float value, int x, int y )
00505 {
00506 array[x][y] = value ;
00507 }
00508
00509 double sinfo_new_doublearray_get_value( double * array, int i )
00510 {
00511 return array[i] ;
00512 }
00513 void sinfo_new_doublearray_set_value( double * array, double value, int i )
00514 {
00515 array[i] = value ;
00516 }
00517
00518 void sinfo_new_destroy_doublearray(double * array)
00519 {
00520 cpl_free( array ) ;
00521 }
00522 double * sinfo_new_doublearray( int size)
00523 {
00524 return (double *) cpl_calloc( size, sizeof (double) ) ;
00525 }
00526
00527 double ** sinfo_new_2Ddoublearray( int size_x, int size_y)
00528 {
00529 double ** retVal ;
00530 int i ;
00531
00532 retVal = (double **) cpl_calloc( size_x, sizeof (double*) ) ;
00533 for ( i = 0 ; i < size_x ; i++ )
00534 {
00535 retVal[i] = (double *) cpl_calloc( size_y, sizeof (double)) ;
00536 }
00537 return retVal ;
00538 }
00539
00540 float ** sinfo_new_2Dfloatarray( int size_x, int size_y)
00541 {
00542 float ** retVal ;
00543 int i ;
00544
00545 retVal = (float **) cpl_calloc( size_x, sizeof (float*) ) ;
00546 for ( i = 0 ; i < size_x ; i++ )
00547 {
00548 retVal[i] = (float *) cpl_calloc( size_y, sizeof (float)) ;
00549 }
00550 return retVal ;
00551 }
00552
00553
00554 int * sinfo_new_intarray( int size)
00555 {
00556 return (int *) cpl_calloc( size, sizeof (int) ) ;
00557 }
00558 void sinfo_new_destroy_intarray(int ** array)
00559 {
00560 cpl_free( *array ) ;
00561 *array=NULL;
00562 }
00563
00564 int sinfo_new_intarray_get_value( int * array, int i )
00565 {
00566 return array[i] ;
00567 }
00568
00569 float sinfo_new_Stats_get_cleanstdev(Stats * stats)
00570 {
00571 return stats -> cleanstdev ;
00572 }
00573 float sinfo_new_Stats_get_cleanmean(Stats * stats)
00574 {
00575 return stats -> cleanmean ;
00576 }
00577
00578 char * sinfo_new_get_basename(const char *filename)
00579 {
00580 char *p ;
00581 p = strrchr (filename, '/');
00582 return p ? p + 1 : (char *) filename;
00583 }
00584
00585
00586
00587 char * sinfo_new_get_rootname(const char * filename)
00588 {
00589 static char path[MAX_NAME_SIZE+1];
00590 char * lastdot ;
00591
00592 if (strlen(filename)>MAX_NAME_SIZE) return NULL ;
00593 memset(path, MAX_NAME_SIZE, 0);
00594 strcpy(path, filename);
00595 lastdot = strrchr(path, '.');
00596 if (lastdot == NULL) return path ;
00597 if ((!strcmp(lastdot, ".fits")) || (!strcmp(lastdot, ".FITS")) ||
00598 (!strcmp(lastdot, ".paf")) || (!strcmp(lastdot, ".PAF")) ||
00599 (!strcmp(lastdot, ".dat")) || (!strcmp(lastdot, ".DAT")) ||
00600 (!strcmp(lastdot, ".fits")) || (!strcmp(lastdot, ".TFITS")) ||
00601 (!strcmp(lastdot, ".ascii")) || (!strcmp(lastdot, ".ASCII")))
00602 {
00603 lastdot[0] = (char)0;
00604 }
00605 return path ;
00606 }
00607
00608
00609
00610
00618
00619 cpl_imagelist * sinfo_new_frameset_to_iset(cpl_frameset * fset)
00620 {
00621 cpl_imagelist * iset=NULL ;
00622 char ** filenames ;
00623 int nfiles=0 ;
00624
00625
00626 if (fset == NULL) return NULL ;
00627
00628
00629 if ((filenames = sinfo_new_frameset_to_filenames(fset, &nfiles)) == NULL) {
00630 sinfo_msg_error( "Cannot get the files names") ;
00631 return NULL ;
00632 }
00633
00634 if ((iset = sinfo_new_imagelist_load_frameset(fset,
00635 CPL_TYPE_FLOAT, 0, 0)) == NULL) {
00636 sinfo_msg_error( "Cannot load *** the image set") ;
00637 sinfo_msg_error((char* ) cpl_error_get_message());
00638
00639 cpl_free(filenames) ;
00640 return NULL ;
00641 }
00642
00643
00644 cpl_free(filenames) ;
00645 return iset ;
00646 }
00647 #include "cpl_imagelist_io.h"
00648
00649
00660
00661 cpl_imagelist *
00662 sinfo_new_imagelist_load_frameset(const cpl_frameset * frameset,
00663 cpl_type type,
00664 int pnum,
00665 int extnum)
00666 {
00667 cpl_image * image = NULL;
00668 cpl_imagelist * imagelist = NULL;
00669 cpl_frame * frame = cpl_frameset_get_first(frameset);
00670 const int nz = cpl_frameset_get_size(frameset);
00671 int i;
00672
00673
00674 cpl_ensure(nz > 0, CPL_ERROR_DATA_NOT_FOUND, NULL);
00675
00676 for (i = 0; frame != NULL; i++, frame = cpl_frameset_get_next(frameset)) {
00677
00678 const char * name = cpl_frame_get_filename(frame);
00679 if (name == NULL) break;
00680
00681
00682 image = cpl_image_load(name, type, pnum, extnum);
00683
00684 if (image == NULL) break;
00685
00686 if (i == 0) {
00687 const int nx = cpl_image_get_size_x(image);
00688 const int ny = cpl_image_get_size_y(image);
00689
00690 if (nx < 1 || ny < 1) break;
00691 imagelist = cpl_imagelist_new();
00692 if (imagelist == NULL) break;
00693 }
00694
00695 if (cpl_imagelist_set(imagelist, image, i)) break;
00696 image = NULL;
00697
00698 }
00699
00700 if (i != nz) {
00701
00702 cpl_image_delete(image);
00703 cpl_imagelist_delete(imagelist);
00704 imagelist = NULL;
00705 }
00706 return imagelist;
00707
00708 }
00709
00720 char ** sinfo_new_frameset_to_filenames(cpl_frameset *set, int *nfiles)
00721 {
00722 char **filenames=NULL;
00723
00724 int nbframes=0;
00725 int i=0;
00726
00727 cpl_frame *curr_frame;
00728
00729 if (set == NULL) {
00730 return NULL;
00731 }
00732
00733 nbframes = cpl_frameset_get_size(set);
00734
00735 if (nbframes < 1) {
00736 return NULL;
00737 }
00738
00739
00740
00741 filenames = cpl_malloc(nbframes * sizeof(char *));
00742
00743 curr_frame = cpl_frameset_get_first(set);
00744 for (i = 0; i < nbframes; i++) {
00745 filenames[i]=(char*) cpl_frame_get_filename(curr_frame);
00746 curr_frame = cpl_frameset_get_next(set);
00747 }
00748
00749
00750
00751
00752
00753 *nfiles = nbframes;
00754
00755 return filenames;
00756
00757 }
00758
00759