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 #ifdef HAVE_CONFIG_H
00028 #include <config.h>
00029 #endif
00030
00031
00036
00037
00038
00042
00043
00044
00045
00046 #include <math.h>
00047
00048 #include <xsh_data_dispersol.h>
00049 #include <xsh_utils.h>
00050 #include <xsh_error.h>
00051 #include <xsh_msg.h>
00052 #include <xsh_pfits.h>
00053 #include <xsh_dfs.h>
00054 #include <cpl.h>
00055 #include <xsh_utils_table.h>
00056 #include <xsh_data_spectralformat.h>
00057 #include <xsh_data_order.h>
00058
00059
00060
00061
00062
00063
00064
00076
00077 xsh_dispersol_list* xsh_dispersol_list_new( int size, int degx, int degy,
00078 xsh_instrument *instrument)
00079 {
00080 xsh_dispersol_list *result = NULL;
00081
00082
00083 XSH_ASSURE_NOT_ILLEGAL( size > 0);
00084 XSH_ASSURE_NOT_NULL( instrument);
00085
00086 XSH_CALLOC( result, xsh_dispersol_list, 1);
00087 result->size = size;
00088 result->degx = degx;
00089 result->degy = degy;
00090 check( result->binx = xsh_instrument_get_binx( instrument));
00091 check( result->biny = xsh_instrument_get_biny( instrument));
00092 XSH_CALLOC( result->list, xsh_dispersol, result->size);
00093 XSH_NEW_PROPERTYLIST( result->header);
00094
00095 cleanup:
00096 if ( cpl_error_get_code() != CPL_ERROR_NONE) {
00097 xsh_dispersol_list_free( &result);
00098 }
00099 return result;
00100 }
00101
00102
00103
00104
00113
00114 xsh_dispersol_list* xsh_dispersol_list_load( cpl_frame *frame,
00115 xsh_instrument *instrument)
00116 {
00117 xsh_dispersol_list *result = NULL;
00118 const char* tablename = NULL;
00119 cpl_table* table = NULL;
00120 int rows, cols ;
00121 int degx, degy;
00122 char coefname[20];
00123 int i;
00124
00125
00126 XSH_ASSURE_NOT_ILLEGAL( frame);
00127
00128
00129 check( tablename = cpl_frame_get_filename(frame));
00130
00131 XSH_TABLE_LOAD( table, tablename);
00132 check( rows = cpl_table_get_nrow( table));
00133 check( cols = cpl_table_get_ncol( table));
00134
00135
00136 check( xsh_get_table_value( table, XSH_DISPERSOL_TABLE_COLNAME_DEGX,
00137 CPL_TYPE_INT, 0, °x));
00138
00139 check( xsh_get_table_value( table, XSH_DISPERSOL_TABLE_COLNAME_DEGY,
00140 CPL_TYPE_INT, 0, °y));
00141
00142 check( result = xsh_dispersol_list_new( rows/2, degx, degy, instrument));
00143
00144 for(i=0; i< rows/2; i++){
00145 cpl_polynomial *lambda_poly = NULL;
00146 cpl_polynomial *slit_poly = NULL;
00147 cpl_size pows[2];
00148 int k,l, absorder = 0;
00149 int j = 2*i;
00150
00151 check( xsh_get_table_value( table, XSH_DISPERSOL_TABLE_COLNAME_ORDER,
00152 CPL_TYPE_INT, j, &absorder));
00153 check( lambda_poly = cpl_polynomial_new( 2));
00154 check( slit_poly = cpl_polynomial_new( 2));
00155 for( k=0; k<= degx; k++){
00156 for( l=0; l<= degy; l++){
00157 double coef_lambda = 0.0, coef_slit = 0.0;
00158 sprintf(coefname,"C%d%d",k,l);
00159 check( xsh_get_table_value( table, coefname,CPL_TYPE_DOUBLE, j,
00160 &coef_lambda));
00161 check( xsh_get_table_value( table, coefname,CPL_TYPE_DOUBLE, j+1,
00162 &coef_slit));
00163 pows[0] =k;
00164 pows[1] = l;
00165 check( cpl_polynomial_set_coeff( lambda_poly, pows, coef_lambda));
00166 check( cpl_polynomial_set_coeff( slit_poly, pows, coef_slit));
00167 }
00168 }
00169 check( xsh_dispersol_list_add( result, i, absorder, lambda_poly, slit_poly));
00170 }
00171
00172 cleanup:
00173 xsh_free_table( &table);
00174 if ( cpl_error_get_code() != CPL_ERROR_NONE) {
00175 xsh_dispersol_list_free( &result);
00176 }
00177 return result;
00178 }
00179
00180
00181
00196
00197 void xsh_dispersol_list_add( xsh_dispersol_list *list, int idx,
00198 int absorder, cpl_polynomial *lambda_poly,
00199 cpl_polynomial *slit_poly)
00200 {
00201
00202 XSH_ASSURE_NOT_NULL( list);
00203 XSH_ASSURE_NOT_NULL( lambda_poly);
00204 XSH_ASSURE_NOT_NULL( slit_poly);
00205 XSH_ASSURE_NOT_ILLEGAL( idx >=0 && idx < list->size);
00206
00207 list->list[idx].absorder = absorder;
00208 list->list[idx].lambda_poly = lambda_poly;
00209 list->list[idx].slit_poly = slit_poly;
00210
00211 cleanup:
00212 return;
00213 }
00214
00215
00216
00223
00224 void xsh_dispersol_list_free( xsh_dispersol_list **list)
00225 {
00226 int i = 0;
00227
00228 if ( list && *list){
00229
00230 for( i = 0; i < (*list)->size; i++){
00231 xsh_free_polynomial( &(*list)->list[i].lambda_poly);
00232 xsh_free_polynomial( &(*list)->list[i].slit_poly);
00233 }
00234 if ((*list)->list){
00235 cpl_free( (*list)->list);
00236 }
00237 xsh_free_propertylist( &((*list)->header));
00238 cpl_free( *list);
00239 *list = NULL;
00240 }
00241 }
00242
00243
00244
00245
00260
00261 cpl_frame* xsh_dispersol_list_to_wavemap( xsh_dispersol_list *list,
00262 cpl_frame *order_frame,
00263 xsh_pre *pre,
00264 xsh_instrument *instr,
00265 const char* tag)
00266
00267 {
00268 cpl_frame *result = NULL;
00269 cpl_image *image = NULL;
00270 xsh_order_list *order_list = NULL;
00271 int nx, ny;
00272 cpl_vector* pos = NULL;
00273 double *data = NULL;
00274 int i, y;
00275 char filename[80];
00276 cpl_propertylist *wavemap_header = NULL;
00277
00278
00279 XSH_ASSURE_NOT_NULL( list);
00280 XSH_ASSURE_NOT_NULL( order_frame);
00281 XSH_ASSURE_NOT_NULL( pre);
00282 XSH_ASSURE_NOT_NULL( instr);
00283
00284
00285 check( order_list = xsh_order_list_load( order_frame, instr));
00286 nx = pre->nx;
00287 ny = pre->ny;
00288 pos = cpl_vector_new(2);
00289
00290
00291 check( image = cpl_image_new( nx, ny, CPL_TYPE_DOUBLE));
00292 check( data = cpl_image_get_data_double( image));
00293 check( wavemap_header = cpl_propertylist_new());
00294
00295
00296 for( i=0; i < list->size; i++) {
00297 int starty, endy;
00298 double lambda_min, lambda_max;
00299 double lambda_starty_lo, lambda_starty_up, lambda_endy_lo, lambda_endy_up;
00300 double dxmin, dxmax ;
00301 int xmin, xmax;
00302 int absorder;
00303
00304
00305 absorder = order_list->list[i].absorder;
00306 check( starty = xsh_order_list_get_starty( order_list, i));
00307
00308 check( endy = xsh_order_list_get_endy( order_list, i));
00309
00310
00311 for ( y=starty ; y<=endy && y<ny; y++){
00312 int x;
00313
00314
00315 check( dxmin = xsh_order_list_eval( order_list,
00316 order_list->list[i].edglopoly,
00317 (double)y ) ) ;
00318 check( dxmax = xsh_order_list_eval( order_list,
00319 order_list->list[i].edguppoly,
00320 (double)y ) ) ;
00321 xmin = ceil( dxmin);
00322 xmax = floor( dxmax);
00323 for ( x=xmin ; x<=xmax && x<nx; x++){
00324 double lambda;
00325
00326 cpl_vector_set(pos, 0, (double)x);
00327 cpl_vector_set(pos, 1, (double)y);
00328
00329
00330 check( lambda = xsh_dispersol_list_eval(
00331 list, list->list[i].lambda_poly, pos));
00332
00333 data[x-1+(y-1)*nx] = (float)lambda;
00334
00335 }
00336
00337 }
00338
00339
00340 check( dxmin = xsh_order_list_eval( order_list,
00341 order_list->list[i].edglopoly,
00342 (double)starty));
00343 check( dxmax = xsh_order_list_eval( order_list,
00344 order_list->list[i].edguppoly,
00345 (double)starty));
00346 xmin = ceil( dxmin);
00347 xmax = floor( dxmax);
00348
00349 cpl_vector_set(pos, 0, (double)xmin);
00350 cpl_vector_set(pos, 1, (double)starty);
00351 check( lambda_starty_lo = xsh_dispersol_list_eval(
00352 list, list->list[i].lambda_poly, pos));
00353
00354 cpl_vector_set(pos, 0, (double)xmax);
00355 check( lambda_starty_up = xsh_dispersol_list_eval(
00356 list, list->list[i].lambda_poly, pos));
00357
00358 check( dxmin = xsh_order_list_eval( order_list,
00359 order_list->list[i].edglopoly,
00360 (double)endy));
00361 check( dxmax = xsh_order_list_eval( order_list,
00362 order_list->list[i].edguppoly,
00363 (double)endy));
00364 xmin = ceil( dxmin);
00365 xmax = floor( dxmax);
00366
00367 cpl_vector_set(pos, 0, (double)xmin);
00368 cpl_vector_set(pos, 1, (double)endy);
00369 check( lambda_endy_lo = xsh_dispersol_list_eval(
00370 list, list->list[i].lambda_poly, pos));
00371
00372 cpl_vector_set(pos, 0, (double)xmax);
00373 check( lambda_endy_up = xsh_dispersol_list_eval(
00374 list, list->list[i].lambda_poly, pos));
00375
00376 if ( lambda_starty_lo < lambda_endy_lo
00377 && lambda_starty_up < lambda_endy_up){
00378 if ( lambda_starty_lo > lambda_starty_up){
00379 lambda_min = lambda_starty_lo;
00380 }
00381 else{
00382 lambda_min = lambda_starty_up;
00383 }
00384 if ( lambda_endy_lo > lambda_endy_up){
00385 lambda_max = lambda_endy_up;
00386 }
00387 else{
00388 lambda_max = lambda_endy_lo;
00389 }
00390 }
00391 else if ( lambda_starty_lo > lambda_endy_lo
00392 && lambda_starty_up > lambda_endy_up){
00393 if ( lambda_starty_lo > lambda_starty_up){
00394 lambda_max = lambda_starty_up;
00395 }
00396 else{
00397 lambda_max = lambda_starty_lo;
00398 }
00399 if ( lambda_endy_lo > lambda_endy_up){
00400 lambda_min = lambda_endy_lo;
00401 }
00402 else{
00403 lambda_min = lambda_endy_up;
00404 }
00405 }
00406 else{
00407 xsh_msg_warning(
00408 "wavelength variation differs between upper(%f %f) and lower edges( %f %f)",
00409 lambda_starty_up, lambda_endy_up, lambda_endy_up, lambda_endy_lo);
00410 lambda_min = 0.0;
00411 lambda_max = 0.0;
00412 }
00413
00414
00415 check( xsh_pfits_set_wavemap_order_lambda_min(
00416 wavemap_header, absorder, lambda_min));
00417 check( xsh_pfits_set_wavemap_order_lambda_max(
00418 wavemap_header, absorder, lambda_max));
00419 }
00420
00421 sprintf(filename,"%s.fits",tag);
00422 check( xsh_pfits_set_pcatg( wavemap_header, tag));
00423
00424 check( cpl_image_save( image, filename, CPL_BPP_IEEE_FLOAT,
00425 wavemap_header, CPL_IO_DEFAULT));
00426
00427 check( result = xsh_frame_product(filename, tag,
00428 CPL_FRAME_TYPE_IMAGE, CPL_FRAME_GROUP_PRODUCT,
00429 CPL_FRAME_LEVEL_TEMPORARY));
00430
00431 cleanup:
00432 if (cpl_error_get_code() != CPL_ERROR_NONE){
00433 xsh_free_frame( &result);
00434 }
00435 xsh_free_image ( &image);
00436 xsh_free_propertylist( &wavemap_header);
00437 xsh_order_list_free( &order_list);
00438 xsh_free_vector( &pos);
00439
00440
00441 return result;
00442 }
00443
00444
00445
00446
00461
00462 cpl_frame* xsh_dispersol_list_to_slitmap( xsh_dispersol_list *list,
00463 cpl_frame *order_frame,
00464 xsh_pre *pre,
00465 xsh_instrument *instr,
00466 const char* tag)
00467 {
00468 cpl_frame *result = NULL;
00469 cpl_image *image = NULL;
00470 cpl_image *ifu_map = NULL;
00471 xsh_order_list *order_list = NULL;
00472 int nx, ny;
00473 cpl_vector* pos = NULL;
00474 double *data = NULL;
00475 int *imap = NULL;
00476
00477 int i, y;
00478 char filename[80];
00479 cpl_vector *med_low_vect = NULL, *med_up_vect = NULL, *med_cen_vect = NULL;
00480 cpl_vector *med_sliclo_vect = NULL, *med_slicup_vect = NULL;
00481 int nb_orders;
00482 double slit_up, slit_low, slit_cen;
00483 cpl_propertylist *slitmap_header = NULL;
00484
00485
00486 XSH_ASSURE_NOT_NULL( list);
00487 XSH_ASSURE_NOT_NULL( order_frame);
00488 XSH_ASSURE_NOT_NULL( pre);
00489 XSH_ASSURE_NOT_NULL( instr);
00490
00491
00492 check( order_list = xsh_order_list_load( order_frame, instr));
00493 nx = pre->nx;
00494 ny = pre->ny;
00495
00496
00497 XSH_ASSURE_NOT_ILLEGAL( order_list->bin_x == pre->binx);
00498 XSH_ASSURE_NOT_ILLEGAL( order_list->bin_x == list->binx);
00499
00500 pos = cpl_vector_new(2);
00501
00502 check( image = cpl_image_new( nx, ny, CPL_TYPE_DOUBLE));
00503 check( data = cpl_image_get_data_double( image));
00504 check( ifu_map = cpl_image_new( nx, ny, CPL_TYPE_INT));
00505 check( imap = cpl_image_get_data_int( ifu_map));
00506 check( slitmap_header = cpl_propertylist_new());
00507
00508 nb_orders = list->size;
00509
00510 med_low_vect = cpl_vector_new( nb_orders);
00511 med_cen_vect = cpl_vector_new( nb_orders);
00512 med_up_vect = cpl_vector_new( nb_orders);
00513 if ( xsh_instrument_get_mode( instr) == XSH_MODE_IFU){
00514 med_sliclo_vect = cpl_vector_new( nb_orders);
00515 med_slicup_vect = cpl_vector_new( nb_orders);
00516 }
00517
00518
00519 for( i=0; i < nb_orders; i++) {
00520 int starty, endy, vect_size;
00521 int absorder;
00522 cpl_vector *low_vect = NULL, *up_vect = NULL, *cen_vect = NULL;
00523 cpl_vector *sliclo_vect = NULL, *slicup_vect = NULL;
00524
00525 check( starty = xsh_order_list_get_starty( order_list, i));
00526 check( endy = xsh_order_list_get_endy( order_list, i));
00527 absorder = order_list->list[i].absorder;
00528
00529 vect_size = endy-starty+1;
00530
00531 low_vect = cpl_vector_new( vect_size);
00532 cen_vect = cpl_vector_new( vect_size);
00533 up_vect = cpl_vector_new( vect_size);
00534
00535 if ( xsh_instrument_get_mode( instr) == XSH_MODE_IFU){
00536 sliclo_vect = cpl_vector_new( vect_size);
00537 slicup_vect = cpl_vector_new( vect_size);
00538 }
00539
00540 for ( y=starty ; y<=endy && y<ny; y++){
00541 double dxmin, dxmax, xcen;
00542 int xmin, xmax, x;
00543 double slit;
00544 double ypos;
00545
00546
00547 check( dxmin = xsh_order_list_eval( order_list,
00548 order_list->list[i].edglopoly,
00549 (double)y ) ) ;
00550 check( dxmax = xsh_order_list_eval( order_list,
00551 order_list->list[i].edguppoly,
00552 (double)y ) ) ;
00553 check( xcen = xsh_order_list_eval( order_list,
00554 order_list->list[i].cenpoly,
00555 (double)y ) ) ;
00556
00557 xmin = ceil( dxmin);
00558 xmax = floor( dxmax);
00559
00560 ypos = (double)y;
00561 cpl_vector_set(pos, 1, ypos);
00562
00563 for ( x=xmin ; x<=xmax && x<nx; x++){
00564 double xpos;
00565
00566 xpos = (double)x;
00567 cpl_vector_set(pos, 0, xpos);
00568 cpl_vector_set(pos, 1, ypos);
00569
00570 check( slit = xsh_dispersol_list_eval( list, list->list[i].slit_poly, pos));
00571
00572
00573 data[x-1+(y-1)*nx] = (float)slit;
00574 imap[x-1+(y-1)*nx] = 1;
00575 }
00576
00577 cpl_vector_set(pos, 0, xmin);
00578 cpl_vector_set(pos, 1, ypos);
00579 check( slit = xsh_dispersol_list_eval( list, list->list[i].slit_poly, pos));
00580 cpl_vector_set( low_vect, y-starty, slit);
00581
00582 cpl_vector_set(pos, 0, xmax);
00583 cpl_vector_set(pos, 1, ypos);
00584 check( slit = xsh_dispersol_list_eval( list, list->list[i].slit_poly, pos));
00585 cpl_vector_set( up_vect, y-starty, slit);
00586
00587 cpl_vector_set(pos, 0, xcen);
00588 cpl_vector_set(pos, 1, ypos);
00589 check( slit = xsh_dispersol_list_eval( list, list->list[i].slit_poly, pos));
00590 cpl_vector_set( cen_vect, y-starty, slit);
00591
00592 if ( xsh_instrument_get_mode( instr) == XSH_MODE_IFU){
00593 xmin = xsh_order_list_eval( order_list,
00594 order_list->list[i].sliclopoly,
00595 (double)y );
00596 xmax = xsh_order_list_eval( order_list,
00597 order_list->list[i].slicuppoly,
00598 (double)y );
00599 cpl_vector_set(pos, 0, xmin);
00600 cpl_vector_set(pos, 1, ypos);
00601 check( slit = xsh_dispersol_list_eval( list, list->list[i].slit_poly, pos));
00602 cpl_vector_set( sliclo_vect, y-starty, slit);
00603
00604 cpl_vector_set(pos, 0, xmax);
00605 cpl_vector_set(pos, 1, ypos);
00606 check( slit = xsh_dispersol_list_eval( list, list->list[i].slit_poly, pos));
00607 cpl_vector_set( slicup_vect, y-starty, slit);
00608 }
00609
00610 }
00611
00612 check( slit_up = cpl_vector_get_median( up_vect));
00613 check( slit_low = cpl_vector_get_median( low_vect));
00614 check( slit_cen = cpl_vector_get_median( cen_vect));
00615
00616
00617 check( xsh_pfits_set_slitmap_order_cen( slitmap_header, absorder, slit_cen));
00618 check( cpl_vector_set( med_cen_vect, i, slit_cen));
00619
00620 if ( slit_up > slit_low){
00621 check( xsh_pfits_set_slitmap_order_edgup( slitmap_header, absorder, slit_up));
00622 check( xsh_pfits_set_slitmap_order_edglo( slitmap_header, absorder, slit_low));
00623 check( cpl_vector_set( med_low_vect, i, slit_low));
00624 check( cpl_vector_set( med_up_vect, i, slit_up));
00625 }
00626 else{
00627 check( xsh_pfits_set_slitmap_order_edgup( slitmap_header, absorder, slit_low));
00628 check( xsh_pfits_set_slitmap_order_edglo( slitmap_header, absorder, slit_up));
00629 check( cpl_vector_set( med_low_vect, i, slit_up));
00630 check( cpl_vector_set( med_up_vect, i, slit_low));
00631 }
00632
00633 if ( xsh_instrument_get_mode( instr) == XSH_MODE_IFU){
00634 check( slit_up = cpl_vector_get_median( slicup_vect));
00635 check( slit_low = cpl_vector_get_median( sliclo_vect));
00636 if ( slit_up > slit_low){
00637 check( xsh_pfits_set_slitmap_order_slicup( slitmap_header, absorder, slit_up));
00638 check( xsh_pfits_set_slitmap_order_sliclo( slitmap_header, absorder, slit_low));
00639 check( cpl_vector_set( med_sliclo_vect, i, slit_low));
00640 check( cpl_vector_set( med_slicup_vect, i, slit_up));
00641 }
00642 else{
00643 check( xsh_pfits_set_slitmap_order_slicup( slitmap_header, absorder, slit_low));
00644 check( xsh_pfits_set_slitmap_order_sliclo( slitmap_header, absorder, slit_up));
00645 check( cpl_vector_set( med_sliclo_vect, i, slit_up));
00646 check( cpl_vector_set( med_slicup_vect, i, slit_low));
00647 }
00648 xsh_free_vector( &slicup_vect);
00649 xsh_free_vector( &sliclo_vect);
00650 }
00651 xsh_free_vector( &up_vect);
00652 xsh_free_vector( &low_vect);
00653 xsh_free_vector( &cen_vect);
00654 }
00655
00656 check( slit_cen = cpl_vector_get_median( med_cen_vect));
00657 check( slit_up = cpl_vector_get_median( med_up_vect));
00658 check( slit_low = cpl_vector_get_median( med_low_vect));
00659
00660 check( xsh_pfits_set_slitmap_median_cen( slitmap_header, slit_cen));
00661 check( xsh_pfits_set_slitmap_median_edgup( slitmap_header, slit_up));
00662 check( xsh_pfits_set_slitmap_median_edglo( slitmap_header, slit_low));
00663
00664 if ( xsh_instrument_get_mode( instr) == XSH_MODE_IFU){
00665 check( slit_up = cpl_vector_get_median( med_slicup_vect));
00666 check( slit_low = cpl_vector_get_median( med_sliclo_vect));
00667 check( xsh_pfits_set_slitmap_median_slicup( slitmap_header, slit_up));
00668 check( xsh_pfits_set_slitmap_median_sliclo( slitmap_header, slit_low));
00669 }
00670
00671 sprintf(filename,"%s.fits",tag);
00672 xsh_pfits_set_pcatg( slitmap_header, tag);
00673 check( cpl_image_save( image, filename, CPL_BPP_IEEE_FLOAT,
00674 slitmap_header, CPL_IO_DEFAULT));
00675
00676 if ( xsh_instrument_get_mode( instr) == XSH_MODE_IFU){
00677 check( cpl_image_save( ifu_map, filename, CPL_BPP_32_SIGNED,
00678 NULL, CPL_IO_EXTEND));
00679 }
00680 check( result = xsh_frame_product(filename,tag,
00681 CPL_FRAME_TYPE_IMAGE, CPL_FRAME_GROUP_PRODUCT,
00682 CPL_FRAME_LEVEL_TEMPORARY));
00683
00684 cleanup:
00685 if (cpl_error_get_code() != CPL_ERROR_NONE){
00686 xsh_free_frame( &result);
00687 }
00688
00689 xsh_free_image ( &image);
00690 xsh_free_image ( &ifu_map);
00691 xsh_free_propertylist( &slitmap_header);
00692 xsh_order_list_free( &order_list);
00693 xsh_free_vector( &pos);
00694 xsh_free_vector( &med_up_vect);
00695 xsh_free_vector( &med_low_vect);
00696 xsh_free_vector( &med_cen_vect);
00697 xsh_free_vector( &med_slicup_vect);
00698 xsh_free_vector( &med_sliclo_vect);
00699
00700 return result;
00701 }
00702
00703
00704
00717
00718 double
00719 xsh_dispersol_list_eval( xsh_dispersol_list *list,
00720 cpl_polynomial *poly,
00721 cpl_vector *pos)
00722 {
00723 double result=0;
00724 double y_bin, y_no_bin;
00725 double x_bin, x_no_bin;
00726
00727 XSH_ASSURE_NOT_NULL( list);
00728 XSH_ASSURE_NOT_NULL( pos);
00729 XSH_ASSURE_NOT_NULL( poly);
00730
00731 check( x_bin = cpl_vector_get(pos, 0));
00732 check( y_bin = cpl_vector_get(pos, 1));
00733
00734 x_no_bin = convert_bin_to_data( x_bin, list->binx);
00735 y_no_bin = convert_bin_to_data( y_bin, list->biny);
00736
00737 check( cpl_vector_set(pos, 0, x_no_bin));
00738 check( cpl_vector_set(pos, 1, y_no_bin));
00739
00740 check( result = cpl_polynomial_eval( poly, pos));
00741
00742 cleanup:
00743 return result;
00744 }
00745
00746
00747
00755
00756 cpl_frame*
00757 xsh_dispersol_list_save( xsh_dispersol_list *list, const char* tag)
00758 {
00759 cpl_frame *result = NULL;
00760 cpl_table *table = NULL;
00761 char coefname[20];
00762 int i,j,k;
00763 int nbcoefs =0;
00764 cpl_size pows[2];
00765 char filename[80];
00766
00767
00768 XSH_ASSURE_NOT_NULL( list);
00769
00770 nbcoefs = (list->degx+1)*(list->degy+1);
00771
00772 check( table = cpl_table_new( XSH_DISPERSOL_TABLE_NBCOL+nbcoefs));
00773 check(cpl_table_set_size(table, list->size*2));
00774
00775 check(
00776 cpl_table_new_column( table, XSH_DISPERSOL_TABLE_COLNAME_AXIS,
00777 CPL_TYPE_STRING));
00778 check(
00779 cpl_table_new_column( table, XSH_DISPERSOL_TABLE_COLNAME_ORDER,
00780 CPL_TYPE_INT));
00781 check(
00782 cpl_table_new_column( table, XSH_DISPERSOL_TABLE_COLNAME_DEGX,
00783 CPL_TYPE_INT));
00784 check(
00785 cpl_table_new_column( table, XSH_DISPERSOL_TABLE_COLNAME_DEGY,
00786 CPL_TYPE_INT));
00787
00788 for(i=0; i<= list->degx; i++){
00789 for(j = 0; j<= list->degy; j++){
00790 sprintf(coefname,"C%d%d",i,j);
00791 check( cpl_table_new_column( table, coefname, CPL_TYPE_DOUBLE));
00792 }
00793 }
00794
00795 for(i=0; i< list->size; i++){
00796 int i2 = i*2;
00797 check(cpl_table_set_string(table,XSH_DISPERSOL_TABLE_COLNAME_AXIS,
00798 i2, XSH_DISPERSOL_AXIS_LAMBDA));
00799 check(cpl_table_set_string(table,XSH_DISPERSOL_TABLE_COLNAME_AXIS,
00800 i2+1, XSH_DISPERSOL_AXIS_SLIT));
00801 check(cpl_table_set_int(table,XSH_DISPERSOL_TABLE_COLNAME_ORDER,
00802 i2, list->list[i].absorder));
00803 check(cpl_table_set_int(table,XSH_DISPERSOL_TABLE_COLNAME_ORDER,
00804 i2+1, list->list[i].absorder));
00805 check(cpl_table_set_int(table,XSH_DISPERSOL_TABLE_COLNAME_DEGX,
00806 i2, list->degx));
00807 check(cpl_table_set_int(table,XSH_DISPERSOL_TABLE_COLNAME_DEGX,
00808 i2+1, list->degx));
00809 check(cpl_table_set_int(table,XSH_DISPERSOL_TABLE_COLNAME_DEGY,
00810 i2, list->degy));
00811 check(cpl_table_set_int(table,XSH_DISPERSOL_TABLE_COLNAME_DEGY,
00812 i2+1, list->degy));
00813 for( j=0; j<= list->degx; j++){
00814 for( k=0; k<= list->degy; k++){
00815 double coef_lambda = 0.0, coef_slit = 0.0;
00816 sprintf(coefname,"C%d%d",j,k);
00817 pows[0] = j;
00818 pows[1] = k;
00819 check(coef_lambda =
00820 cpl_polynomial_get_coeff( list->list[i].lambda_poly, pows));
00821 check(cpl_table_set_double( table, coefname, i2, coef_lambda));
00822 check(coef_slit =
00823 cpl_polynomial_get_coeff( list->list[i].slit_poly, pows));
00824 check(cpl_table_set_double( table, coefname, i2+1, coef_slit));
00825 }
00826 }
00827 }
00828
00829 sprintf(filename, "%s.fits",tag);
00830 check( xsh_pfits_set_pcatg(list->header,tag));
00831 check(cpl_table_save(table, list->header, NULL,
00832 filename, CPL_IO_DEFAULT));
00833
00834 check( result= xsh_frame_product( filename, tag,
00835 CPL_FRAME_TYPE_TABLE,
00836 CPL_FRAME_GROUP_PRODUCT,
00837 CPL_FRAME_LEVEL_TEMPORARY));
00838
00839 cleanup:
00840 xsh_free_table( &table);
00841 if ( cpl_error_get_code() != CPL_ERROR_NONE) {
00842 }
00843 return result;
00844 }
00845
00846