diff options
Diffstat (limited to 'funtools/funcopy.c')
-rw-r--r-- | funtools/funcopy.c | 364 |
1 files changed, 364 insertions, 0 deletions
diff --git a/funtools/funcopy.c b/funtools/funcopy.c new file mode 100644 index 0000000..dc2861e --- /dev/null +++ b/funtools/funcopy.c @@ -0,0 +1,364 @@ +/* + * Copyright (c) 1999-2003 Smithsonian Astrophysical Observatory + */ + +#include <funtoolsP.h> + +/* + * + * private routines + * + */ + +/* + * + * semi-public routines -- used by other fun routines, but not by users + * + */ + +#ifdef ANSI_FUNC +void +_FunCopyBinDelete(FITSHead header, int n) +#else +void _FunCopyBinDelete(header, n) + FITSHead header; + int n; +#endif +{ + int i; + char *s; + FITSCard card; + + /* delete unwanted fields */ + ft_headdel(header, "TFIELDS", 0); + /* delete checksum header params */ + if( (s=ft_headgets(header, "CHECKSUM", 0, NULL, &card)) && card){ + ft_headdel(header, "CHECKSUM", 0); + xfree(s); + } + if( (s=ft_headgets(header, "DATASUM", 0, NULL, &card)) && card){ + ft_headdel(header, "DATASUM", 0); + xfree(s); + } + for(i=1; i<=n; i++){ + if( (s=ft_headgets(header, "TFORM", i, NULL, &card)) && card){ + ft_headdel(header, "TFORM", i); + xfree(s); + } + if( (s=ft_headgets(header, "TTYPE", i, NULL, &card)) && card){ + ft_headdel(header, "TTYPE", i); + xfree(s); + } + if( (s=ft_headgets(header, "TUNIT", i, NULL, &card)) && card){ + ft_headdel(header, "TUNIT", i); + xfree(s); + } + if( (s=ft_headgets(header, "TDISP", i, NULL, &card)) && card){ + ft_headdel(header, "TDISP", i); + xfree(s); + } + if( (s=ft_headgets(header, "TSCAL", i, NULL, &card)) && card ){ + ft_headdel(header, "TSCAL", i); + xfree(s); + } + if( (s=ft_headgets(header, "TZERO", i, NULL, &card)) && card ){ + ft_headdel(header, "TZERO", i); + xfree(s); + } + if( (s=ft_headgets(header, "TNULL", i, NULL, &card)) && card ){ + ft_headdel(header, "TNULL", i); + xfree(s); + } + if( (s=ft_headgets(header, "TLMIN", i, NULL, &card)) && card ){ + ft_headdel(header, "TLMIN", i); + xfree(s); + } + if( (s=ft_headgets(header, "TLMAX", i, NULL, &card)) && card ){ + ft_headdel(header, "TLMAX", i); + xfree(s); + } + if( (s=ft_headgets(header, "TDMIN", i, NULL, &card)) && card ){ + ft_headdel(header, "TDMIN", i); + xfree(s); + } + if( (s=ft_headgets(header, "TDMAX", i, NULL, &card)) && card ){ + ft_headdel(header, "TDMAX", i); + xfree(s); + } + if( (s=ft_headgets(header, "TDBIN", i, NULL, &card)) && card ){ + ft_headdel(header, "TDBIN", i); + xfree(s); + } + if( (s=ft_headgets(header, "TCTYP", i, NULL, &card)) && card ){ + ft_headdel(header, "TCTYP", i); + xfree(s); + } + if( (s=ft_headgets(header, "TCRVL", i, NULL, &card)) && card ){ + ft_headdel(header, "TCRVL", i); + xfree(s); + } + if( (s=ft_headgets(header, "TCDLT", i, NULL, &card)) && card ){ + ft_headdel(header, "TCDLT", i); + xfree(s); + } + if( (s=ft_headgets(header, "TCRPX", i, NULL, &card)) && card ){ + ft_headdel(header, "TCRPX", i); + xfree(s); + } + if( (s=ft_headgets(header, "TCROT", i, NULL, &card)) && card ){ + ft_headdel(header, "TCROT", i); + xfree(s); + } + } +} + +/* + * + * _FunCopy2ImageHeader -- copy FITS header, transforming it to primary image. + * This beastly routine is used by FunImagePut to make ordinary life easier. + * + */ +#ifdef ANSI_FUNC +int +_FunCopy2ImageHeader(Fun from, Fun to) +#else +int _FunCopy2ImageHeader(from, to) + Fun from; /* input FITS header */ + Fun to; /* output FITS header */ +#endif +{ + int simple=1; + int i, j; + int cols; + char *s; + char tbuf[SZ_LINE]; + double dval; + double crpix1, crpix2; + double cdelt1, cdelt2; + double ltv1, ltv2; + double ltm1_1, ltm2_1, ltm1_2, ltm2_2; + FITSCard card; + + /* make sure we have something to work with */ + if( !from ) + return 0; + + /* make a copy of the header */ + to->header = ft_headcopy(from->header); + + /* make sure this is a primary FITS file */ + ft_cardfmt((FITSCard)to->header->cards, + "SIMPLE", 0, FT_LOGICAL, &simple, 0, "Standard FITS"); + + /* add/edit image-related keywords based on the orig binning values */ + ft_headseti(to->header, + "NAXIS", 0, from->odims, "number of axes", 1); + ft_headseti(to->header, + "NAXIS", 1, from->odim1, "x axis dimension", 1); + if( from->odims >= 2 ){ + ft_headseti(to->header, + "NAXIS", 2, from->odim2, "y axis dimension", 1); + + } + else{ + ft_headdel(to->header, "NAXIS", 2); + } + ft_headseti(to->header, + "BITPIX", 0, from->obitpix, "bits/pixel", 1); + /* synchronize the header and the cards after any changes */ + ft_syncdata(to->header); + + /* get dimensions of the new image back into the struct */ + to->dim1 = ft_headgeti(to->header, "NAXIS", 1, 0, &card); + to->dim2 = ft_headgeti(to->header, "NAXIS", 2, 0, &card); + /* determine dimensionality */ + if( card ) + to->dims = 2; + else + to->dims = 1; + /* get bitpix */ + to->bitpix = ft_headgeti(to->header, "BITPIX", 0, 0, &card); + + /* now get image size based on these params */ + _FunImageSize(to); + + /* delete checksum header params */ + if( ft_headgets(to->header, "CHECKSUM", 0, NULL, &card) && card) + ft_headdel(to->header, "CHECKSUM", 0); + if( ft_headgets(to->header, "DATASUM", 0, NULL, &card) && card) + ft_headdel(to->header, "DATASUM", 0); + + /* now perform the transformation based on input type */ + switch(from->type){ + case FUN_TABLE: + case FUN_EVENTS: + /* delete extraneous bintable keyword */ + /* but leave the XTENSION card so we can overwrite it later */ + /* ft_headdel(to->header, "XTENSION", 0); */ + ft_headdel(to->header, "PCOUNT", 0); + ft_headdel(to->header, "GCOUNT", 0); + if( (s=ft_headgets(to->header, "EXTNAME", 0, NULL, &card)) && card){ + ft_headdel(to->header, "EXTNAME", 0); + xfree(s); + } + if( (ft_headgeti(to->header, "EXTVER", 0, 0, &card)!=0) && card) + ft_headdel(to->header, "EXTVER", 0); + if( (ft_headgeti(to->header, "EXTLEVEL", 0, 0, &card)!=0) && card) + ft_headdel(to->header, "EXTLEVEL", 0); + /* process table columns */ + cols = ft_headgeti(to->header, "TFIELDS", 0, 0, &card); + /* transform WCS cards for the two columns on which we are binning */ + /* this assumes the TFORMnnn are in order ... hmmm! */ + for(j=1; j<=2; j++){ + if( (i=from->bin[j-1]+1) >= 1 ){ + if( (s=ft_headgets(from->header, "TCTYP", i, NULL, &card)) && card ){ + ft_headapps(to->header, "CTYPE", j, s, NULL); + xfree(s); + } + if( (s=ft_headgets(from->header, "TCRVL", i, NULL, &card)) && card ){ + ft_headappv(to->header, "CRVAL", j, s, NULL); + xfree(s); + } + if((dval=ft_headgetr(from->header, "TCDLT", i, 0.0, &card)) && card){ + /* factor in binsize */ + if( from->cols[i-1] && (from->cols[i-1]->binsiz > 0) ) + dval *= from->cols[i-1]->binsiz; + ft_headsetr(to->header, "CDELT", j, dval, 7, NULL, 1); + } + if( (s=ft_headgets(from->header, "TCROT", i, NULL, &card)) && card ){ + ft_headappv(to->header, "CROTA", j, s, NULL); + xfree(s); + } + if((dval=ft_headgetr(from->header, "TCRPX", i, 0.0, &card)) && card){ + /* convert WCS center value from physical to image coords */ + if( from->cols[i-1] ){ + dval = tlp2i(dval, + from->cols[i-1]->tlmin, from->cols[i-1]->binsiz, 'D'); + } + ft_headsetr(to->header, "CRPIX", j, dval, 7, NULL, 1); + } + } + } + /* delete unwanted fields */ + _FunCopyBinDelete(to->header, cols); + break; + case FUN_IMAGE: + /* delete all signs of higher axes */ + for(i=3; i<=ft_naxes(to->header); i++) + ft_headdel(to->header, "NAXIS", i); + /* if bscale and bzero were applied to the input, delete them */ + if( from->scaled & FUN_SCALE_APPLIED ){ + if( ft_headfind(to->header, "BSCALE", 0, 0) ) + ft_headdel(to->header, "BSCALE", 0); + if( ft_headfind(to->header, "BZERO", 0, 0) ) + ft_headdel(to->header, "BZERO", 0); + } + /* if an image section is blocked and we are summing, + we have to take that into account when setting the bzero factor */ + else if( (ft_bzero(to->header) != 0.0) && + (from->block != 1) && (from->btype == FUN_SUM) ){ + dval = ft_bzero(to->header) * from->block * from->block; + ft_headsetr(to->header, "BZERO", 0, dval, 7, "", 1); + } + /* deal with the IRAF DETSEC keyword */ + if( (s=ft_headgets(from->header, "DATASEC", 0, NULL, &card)) && card ){ + int dx0, dx1, dy0, dy1; + char s1[SZ_LINE], s2[SZ_LINE], s3[SZ_LINE], s4[SZ_LINE]; + char *t=s; + if( *t == '[' ) t++; + if(sscanf(t, "%[-0-9.*] : %[-0-9.*] , %[-0-9.*] : %[-0-9.*]", + s1, s2, s3, s4) == 4){ + dx0 = MAX(1,atoi(s1)-from->x0+1); + dx1 = MIN(from->odim1,atoi(s2)-from->x0+1); + dy0 = MAX(1,atoi(s3)-from->y0+1); + dy1 = MIN(from->odim2,atoi(s4)-from->y0+1); + snprintf(s1, SZ_LINE, "[%d:%d,%d:%d]", dx0, dx1, dy0, dy1); + ft_headsets(to->header, "DATASEC", 0, s1, NULL, 1); + } + xfree(s); + } + break; + case FUN_ARRAY: + default: + break; + } + + /* update WCS values */ + /* re-set the tangent plane point */ + crpix1 = ft_headgetr(to->header, "CRPIX", 1, 0.0, &card); + if( card != NULL ){ + crpix1 = (crpix1 + 1.0 - from->x0 - 0.5)/from->block + 0.5; + ft_headsetr(to->header, "CRPIX", 1, crpix1, 7, "reference point", 1); + } + crpix2 = ft_headgetr(to->header, "CRPIX", 2, 0.0, &card); + if( card != NULL ){ + crpix2 = (crpix2 + 1.0 - from->y0 - 0.5)/from->block + 0.5; + ft_headsetr(to->header, "CRPIX", 2, crpix2, 7, "reference point", 1); + } + + /* update degrees/pixel */ + cdelt1 = ft_headgetr(to->header, "CDELT", 1, 0.0, &card); + if( card != NULL ){ + cdelt1 *= from->block; + ft_headsetr(to->header, "CDELT", 1, cdelt1, 7, "degrees/pixel", 1); + } + cdelt2 = ft_headgetr(to->header, "CDELT", 2, 0.0, &card); + if( card != NULL ){ + cdelt2 *= from->block; + ft_headsetr(to->header, "CDELT", 2, cdelt2, 7, "degrees/pixel", 1); + } + + /* WCS matrix-style */ + for(i=1; i<=2; i++){ + for(j=1; j<=2; j++){ + snprintf(tbuf, SZ_LINE, "CD%d_%d", i, j); + dval = ft_headgetr(to->header, tbuf, 0, 0.0, &card); + if( card != NULL ){ + dval *= from->block; + ft_headsetr(to->header, tbuf, 0, dval, 7, "WCS matrix value", 1); + } + } + } + + /* re-set the tangent plane point for IRAF matrix */ + /* for images, ltv values are in file */ + if( from->odims >= 2 ){ + if( from->header->image ){ + ltv1 = ft_headgetr(from->header, "LTV", 1, 0.0, &card); + ltv2 = ft_headgetr(from->header, "LTV", 2, 0.0, &card); + } + /* for tables, ltv is related to tlmin value of binned columns */ + else{ + ltv1 = 1.0 - (int)(from->min1 + 0.5); + ltv2 = 1.0 - (int)(from->min2 + 0.5); + } + ltv1 = (ltv1 + 1.0 - from->x0 - 0.5 )/(from->block * from->binsiz1) + 0.5; + ltv2 = (ltv2 + 1.0 - from->y0 - 0.5 )/(from->block * from->binsiz2) + 0.5; + ft_headsetr(to->header, "LTV", 1, ltv1, 7, "IRAF ref. point", 1); + ft_headsetr(to->header, "LTV", 2, ltv2, 7, "IRAF ref. point", 1); + + /* IRAF matrix values */ + ltm1_1 = ft_headgetr(from->header, "LTM1_1", 0, 1.0, &card); + ltm2_1 = ft_headgetr(from->header, "LTM2_1", 0, 0.0, &card); + ltm1_2 = ft_headgetr(from->header, "LTM1_2", 0, 0.0, &card); + ltm2_2 = ft_headgetr(from->header, "LTM2_2", 0, 1.0, &card); + ltm1_1 /= (from->block * from->binsiz1); + ltm2_1 /= (from->block * from->binsiz2); + ltm1_2 /= (from->block * from->binsiz1); + ltm2_2 /= (from->block * from->binsiz2); + ft_headsetr(to->header, "LTM1_1", 0, ltm1_1, 7, "IRAF matrix value", 1); + ft_headsetr(to->header, "LTM2_1", 0, ltm2_1, 7, "IRAF matrix value", 1); + ft_headsetr(to->header, "LTM1_2", 0, ltm1_2, 7, "IRAF matrix value", 1); + ft_headsetr(to->header, "LTM2_2", 0, ltm2_2, 7, "IRAF matrix value", 1); + } + + return(1); +} + +/* + * + * public routines + * + */ + +/* none */ |