IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 23924


Ignore:
Timestamp:
Apr 18, 2009, 11:11:15 AM (17 years ago)
Author:
eugene
Message:

import updates from camera group (fits table outputs, fixes to persistent trails); modify import to avoid libpsf (fortran)

Location:
trunk/extsrc/gpcsw
Files:
7 added
16 edited

Legend:

Unmodified
Added
Removed
  • trunk/extsrc/gpcsw/Makefile.in

    r23527 r23924  
    1515        svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/ThisIsTopLevel gpcsrc/ThisIsTopLevel
    1616        svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/missing_protos.h gpcsrc/missing_protos.h
    17 # remove existing directories
    18         rm -rf gpcsrc/fits
    19         rm -rf gpcsrc/analysis
    2017# required subdirs
    21         svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/analysis/libpscoords gpcsrc/analysis/libpscoords
    22         svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/analysis/libpsf gpcsrc/analysis/libpsf
    23         svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/fits/libfh gpcsrc/fits/libfh
    24         svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/fits/libfhreg gpcsrc/fits/libfhreg
    25         svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/fits/burntool gpcsrc/fits/burntool
     18        svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/analysis/libpscoords gpcsrc/analysis/libpscoords.new
     19        svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/analysis/libpsf gpcsrc/analysis/libpsf.new
     20        svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/fits/libfh gpcsrc/fits/libfh.new
     21        svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/fits/libfhreg gpcsrc/fits/libfhreg.new
     22        svn export https://svn.ifa.hawaii.edu/gpc/repos/sw/trunk/fits/burntool gpcsrc/fits/burntool.new
     23# replace new files with those from the imported directories
     24        cp -rf gpcsrc/analysis/libpscoords.new/* gpcsrc/analysis/libpscoords
     25        cp -rf gpcsrc/analysis/libpsf.new/* gpcsrc/analysis/libpsf
     26        cp -rf gpcsrc/fits/libfh.new/* gpcsrc/fits/libfh
     27        cp -rf gpcsrc/fits/libfhreg.new/* gpcsrc/fits/libfhreg
     28        cp -rf gpcsrc/fits/burntool.new/* gpcsrc/fits/burntool
     29# remove the temp directory
     30        rm -rf gpcsrc/analysis/libpscoords.new
     31        rm -rf gpcsrc/analysis/libpsf.new
     32        rm -rf gpcsrc/fits/libfh.new
     33        rm -rf gpcsrc/fits/libfhreg.new
     34        rm -rf gpcsrc/fits/burntool.new
    2635# replace the standard gpc/Make.Common with our repaired version (g77 -> gfortran)
    2736        cp -f Make.Common.fixed gpcsrc/Make.Common
  • trunk/extsrc/gpcsw/gpcsrc/Make.Common

    r23490 r23924  
    134134TAR      = gtar
    135135AR       = ar
    136 FC       = gfortran
     136FC       = g77
    137137CC       = gcc
    138138LD       = gcc
     
    673673        @-rm -f $(DIR_BIN)/*.TEXT_FILE_BUSY 2> /dev/null || /bin/true
    674674        @ln $< $@
    675         # @-setuidinst $@
     675        @-setuidinst $@
    676676
    677677$(DIR_BIN)/%-$(VERSION): scripts/%.sh
  • trunk/extsrc/gpcsw/gpcsrc/analysis/libpscoords/pscoords.h

    r23490 r23924  
    1 /* Three-quarter-assed pscoords.h until jt writes me a real one */
     1/* Header file for pscoords.c */
    22#define NINT(x) (x<0?(int)((x)-0.5):(int)((x)+0.5))
    33
     
    128128const PSC_OFFROT_T psc_otaoff[PSC_NX*PSC_NY];
    129129
     130#define PS_scale 38.860 /* Default PS plate scale [um/arcsec] */
     131#define PS_d3  1.49e-10 /* Default PS distortion [arcsec^-2] */
     132#define PS_airdens 0.71 /* Default PS1 qir density (Haleakala) */
     133
    130134/**********************/
    131135/* Normal Application */
     
    147151 */
    148152
    149 
    150153/**************/
    151154/* Prototypes */
  • trunk/extsrc/gpcsw/gpcsrc/analysis/libpsf/psf.c

    r23490 r23924  
    207207/* Load up extras */
    208208   if(psfextra != NULL) {
     209/* Sanity check against no error but crazy fit! */
     210      psfextra->xfw = xfwhm*binx;
     211      psfextra->yfw = yfwhm*biny;
    209212/* Note that these major/minor things are messed up if binx!=biny */
    210213      psfextra->majfw = wpar[9]*binx;
  • trunk/extsrc/gpcsw/gpcsrc/fits/burntool/Makefile

    r23526 r23924  
    55CCDEFS+=$(VERSIONDEFS)
    66CCLINK += -lm
    7 
    87# $(EXECNAME): $(OBJS) libfh.a libpsf.a libpscoords.a
    98$(EXECNAME): $(OBJS) libfh.a libpscoords.a
     
    1312# Dependencies by Make.Common $Revision: 2.17 $
    1413
    15 $(OBJ)/basement.o: basement.c
    16 
    1714$(OBJ)/burnfix.o: burnfix.c burntool.h burnparams.h
    1815
    1916$(OBJ)/burntool.o: burntool.c fh/fh.h fhreg/general.h \
    2017  fhreg/macros.h fhreg/gpc_detector.h \
    21   burntool.h burnparams.h
     18  burntool.h burnparams.h persist_fits.h
    2219
    2320$(OBJ)/burnutils.o: burnutils.c burntool.h burnparams.h
    2421
    2522$(OBJ)/cheapfits.o: cheapfits.c
     23
     24$(OBJ)/persist_fits.o: persist_fits.c fh/fh.h burntool.h \
     25  burnparams.h persist_fits.h
    2626
    2727$(OBJ)/persistfix.o: persistfix.c burntool.h burnparams.h
     
    3434$(OBJ)/psfstats.o: psfstats.c burntool.h psf/psf.h
    3535
     36$(OBJ)/sort.o: sort.c
     37
    3638$(OBJ)/stardetect.o: stardetect.c burntool.h burnparams.h
    3739
  • trunk/extsrc/gpcsw/gpcsrc/fits/burntool/burnparams.h

    r23490 r23924  
    11/* A global buffer into which routines can accumulate detections */
     2#ifndef _INCLUDED_burnparams_
     3#define _INCLUDED_burnparams_
    24
    35/* '#define EXTERN' in just one file, burntool.c, to declare variables */
     
    2830EXTERN int MAX_READ_NOISE;      /* Maximum believable read noise (ADU) */
    2931EXTERN double MIN_EADU;         /* Minimum believable e/ADU */
     32EXTERN int SAT4SURE;            /* Ignore pixels above for noise estimate */
     33EXTERN double MIN_BLAST_PASS;   /* Allow blasted cells if they have BIG satfrac */
     34EXTERN double MAX_BLAST_PASS;   /* But not if it's all wiped out! */
    3035
    3136EXTERN int BURN_THRESH  ;       /* Threshold for onset of burning */
     
    4853
    4954EXTERN double NEGLIGIBLE_TRAIL; /* Don't sweat less than this * sigma */
     55EXTERN int EXPIRE_TRAIL_TIME;   /* Expire trails after this interval */
    5056
    51 
     57#endif /* _INCLUDED_burnparams_ */
  • trunk/extsrc/gpcsw/gpcsrc/fits/burntool/burntool.c

    r23552 r23924  
    2020#define EXTERN  /* Define EXTERN to declare variables in params.h */
    2121#include "burnparams.h"
     22#include "persist_fits.h"
    2223
    2324int
     
    2829   HeaderUnit ehu;
    2930   const char* ifilename = "-";
    30    const char* ofilename = NULL;
    3131   char extname[FH_MAX_STRLEN+1], otaposn[FH_MAX_STRLEN+1];
    3232   int cellmask[MAXCELL];
     
    3434   int nextend, cellxy, cell, cellcode;
    3535   int ext, update, restore, psfsize, psfavg;
    36    int outfd;
    3736   IMTYPE *buf;
    38    const char *burnfile=NULL,  *persistfile=NULL;
     37   const char *burnfile=NULL,  *persistfile=NULL, *persistfitsfile=NULL;
     38   const char *deltablefitsfile=NULL;
    3939   const char *psffile=NULL, *psfstatfile=NULL;
    4040   CELL OTA[MAXCELL];   /* Cell structure for entire OTA */
    4141
    42    if (argc == 1)  {
    43      syntax(argv[0]);
    44      exit(EXIT_SUCCESS);
    45    }
    46 
    4742   if(argc > 1 && strncmp(argv[1], "help", 4) == 0) {
    48       syntax(argv[0]);
    49       exit(EXIT_SUCCESS);
    50    }
    51 
    52    if(argc > 1 && strncmp(argv[1], "--help", 6) == 0) {
    53       syntax(argv[0]);
    54       exit(EXIT_SUCCESS);
    55    }
    56 
    57    if(argc > 1 && strncmp(argv[1], "-h", 6) == 0) {
    5843      syntax(argv[0]);
    5944      exit(EXIT_SUCCESS);
     
    7156      fprintf(stderr, "\rerror: `%s' is not a multi-extension FITS\n",
    7257              ifilename);
     58#ifndef JT2DHACK
    7359      exit(EXIT_FAILURE);
     60#endif
    7461   }
    7562
     
    9178   MAX_READ_NOISE = 20;         /* Maximum believable read noise (ADU) */
    9279   MIN_EADU = 0.3;              /* Minimum believable e/ADU */
     80   SAT4SURE  = 60000;           /* Ignore pixels above for noise estimate */
     81   MIN_BLAST_PASS = 0.1;        /* Allow blasted cells if they have BIG satfrac */
     82   MAX_BLAST_PASS = 0.9;        /* But not if it's all wiped out! */
    9383
    9484   BURN_THRESH  = 30000;        /* Threshold for onset of burning */
     
    10898
    10999   NEGLIGIBLE_TRAIL = 0.5;      /* Don't sweat less than this * sigma */
    110 
     100   EXPIRE_TRAIL_TIME = 2000;    /* Expire a persist after this [sec] */
    111101
    112102/* Parse the args */
     
    155145         burnfile = argv[i] + 9;
    156146
    157 /* Output file for burn streaks */
     147/* Output text file for burn streaks */
    158148      } else if(strncmp(argv[i], "out=", 4) == 0) {     /* out=fname */
    159149         burnfile = argv[i] + 4;
    160150
    161 # if (0)
    162 /* XXX disable for now: is not yet working */
    163 /* Alternate Output file for deburned image */
    164       } else if(strncmp(argv[i], "outimage=", 9) == 0) {        /* outimage=fname */
    165         ofilename = argv[i] + 9;
    166 # endif
    167 
    168 /* Input file for previous burn persistence streaks */
    169       } else if(strncmp(argv[i], "trailin=", 8) == 0) { /* trailin=fname */
     151/* Input text file for previous burn persistence streaks */
     152      } else if(strncmp(argv[i], "infits=", 8) == 0) {  /* in=fname */
    170153         persistfile = argv[i] + 8;
     154
     155/* Same thing, but information is stored in tables in a FITS file. */
     156      } else if(strncmp(argv[i], "trailinfits=", 8) == 0) { /* infits=fname */
     157         persistfitsfile = argv[i] + 12;
    171158
    172159/* Input file for previous burn persistence streaks */
    173160      } else if(strncmp(argv[i], "in=", 3) == 0) {      /* in=fname */
    174161         persistfile = argv[i] + 3;
     162
     163/* Same thing, but information is stored in tables in a FITS file. */
     164      } else if(strncmp(argv[i], "trailinfits=", 8) == 0) { /* trailin=fname */
     165         persistfitsfile = argv[i] + 12;
    175166
    176167/* Output file for PSF gallery */
     
    244235      } else if(strncmp(argv[i], "thrpsf=", 7) == 0) {  /* thrpsf=thresh */
    245236         if(sscanf(argv[i]+7, "%d", &PSF_THRESH) != 1) {
     237            fprintf(stderr, "\rerror: cannot get value from `%s'\n", argv[i]);
     238            exit(EXIT_FAILURE);
     239         }
     240
     241      } else if(strncmp(argv[i], "expire=", 7) == 0) {  /* expire=nsec */
     242         if(sscanf(argv[i]+7, "%d", &EXPIRE_TRAIL_TIME) != 1) {
    246243            fprintf(stderr, "\rerror: cannot get value from `%s'\n", argv[i]);
    247244            exit(EXIT_FAILURE);
     
    295292         }
    296293
     294      } else if(strncmp(argv[i], "deltables=", 8) == 0) { /* trailin=fname */
     295         deltablefitsfile = argv[i] + 10;       
     296
    297297      } else if(strncmp(argv[i], "help", 4) == 0) {     /* help output */
    298298         syntax(argv[0]);
     
    306306   }
    307307
     308   /* If we're told to remove the table from a FITS file, then
     309    * do nothing else. */
     310   if(deltablefitsfile) {
     311     fprintf(stderr,
     312             "\rRemoving burn tables from %s, all other options ignored.\n",
     313             ifilename);
     314     if(persist_fits_remove_tables(ihu, deltablefitsfile) != FH_SUCCESS) {
     315       exit(EXIT_FAILURE);
     316     }
     317     else {
     318       exit(EXIT_SUCCESS);
     319     }
     320   }
     321   
     322   /* If there is no other persistence info supplied, try getting
     323    * it from the input FITS file. */
    308324   if(restore && persistfile == NULL) {
    309       fprintf(stderr, "\rerror: must specify an input file for restore\n");
    310       exit(EXIT_FAILURE);
     325     if(persistfitsfile == NULL) persistfitsfile = ifilename;
    311326   }
    312327
     
    322337      if(persist_read(OTA, persistfile)) exit(EXIT_FAILURE);
    323338   }
     339   else if(persistfitsfile != NULL) {
     340     if(persist_fits_read(OTA, persistfitsfile) != FH_SUCCESS) exit(EXIT_FAILURE);
     341   }
    324342
    325343/* Which OTA is this??? */
     
    333351   }
    334352
    335    if (ofilename) {
    336      outfd = creat (ofilename, S_IRUSR | S_IWUSR | S_IRGRP | S_IWGRP | S_IROTH | S_IWOTH);
    337      if (outfd == -1) {
    338          fprintf(stderr, "\rerror: Cannot open %s for output\n", ofilename);
    339          exit(EXIT_FAILURE);
    340       }
    341 
    342      if (fh_write(ihu, outfd)) {
    343          fprintf(stderr, "\rerror: Trouble writing PHU to %s\n", ofilename);
    344          exit(EXIT_FAILURE);
    345       }
     353   /* Warn that what you're doing might not be a good idea. */
     354   {
     355     fh_bool burn_applied;
     356     
     357     if(fh_get_bool(ihu, PHU_NAME_BURN_APPLIED,
     358                    &burn_applied) != FH_SUCCESS) {
     359        if(VERBOSE > 0) {
     360           fprintf(stderr,
     361                   "warning: Unable to determine whether burn correction "
     362                   "already applied - unable to find %s in primary header\n",
     363                   PHU_NAME_BURN_APPLIED);
     364        }
     365        burn_applied = FH_FALSE;
     366     }
     367     if(restore && (burn_applied == FH_FALSE)) {
     368       fprintf(stderr,
     369               "warning: Restoring old burns, but header indicates no burns previously corrected.\n");             
     370     }
     371     else if (update && (burn_applied == FH_TRUE)) {
     372       fprintf(stderr,
     373               "warning: Applying burn correction, but header indicates burns previously corrected.\n");             
     374       }
    346375   }
    347376
    348377/* Look at all the MEF's extensions */
     378#ifndef JT2DHACK
    349379   for (ext = 1; ext <= nextend; ext++) {
     380#else
     381   for (ext = 1; ext <= MAX(1,nextend); ext++) {
     382#endif
    350383      int naxis, naxis1, naxis2, naxis3;
    351384      int prescan1, ovrscan1, ovrscan2, pontime;
     
    353386      char xtension[FH_MAX_STRLEN + 1];
    354387
     388#ifdef JT2DHACK
     389      if(nextend == 0) {
     390         ehu = ihu;
     391         sprintf(extname, "%s", "xy00");
     392      } else {
     393#endif
    355394      if (!(ehu = fh_ehu(ihu, ext)) ||
    356395          fh_get_str(ehu, "EXTNAME", extname, sizeof(extname)) != FH_SUCCESS) {
     
    377416         fprintf(stderr,
    378417                 "warning: Skipping non-image extension %s\n", extname);
     418#ifdef JT2DHACK
     419      }
     420#endif
    379421/*
    380422 * %%% When reading from a pipe, does the data still need
     
    448490      if (fh_get_BZERO(ehu, &bzero_d) == FH_SUCCESS) BZERO = bzero_d;
    449491
     492      if(VERBOSE & VERB_NORM) {
     493         printf("nx=%d ny=%d prex=%d postx=%d posty=%d BZERO=%.1f\n",
     494                naxis1, naxis2, prescan1, ovrscan1, ovrscan2, BZERO);
     495      }
     496
    450497      buf = (IMTYPE*)malloc(naxis1*naxis2*naxis3*sizeof(short));
    451498      if (fh_read_padded_image(ehu, fh_file_desc(ehu), buf,
     
    471518/* Get bias and sky levels */
    472519         err = cell_stats(nx, ny, naxis1, naxis2, imbuf, OTA+cell);
     520         if(err) {
     521            if(VERBOSE > 0) {
     522               fprintf(stderr, "logonly: error getting bias/sky/rms for cell %d, skipping...\n", cell);
     523            }
     524            continue;
     525         }
     526
     527/* Does this cell get a pass because it's heavily blasted? */
     528/* Turn it in to a relaxation factor for the "read noise" */
     529         i = OTA[cell].satfrac > MIN_BLAST_PASS &&
     530             OTA[cell].satfrac < MAX_BLAST_PASS &&
     531             OTA[cell].sky < TRAIL_THRESH ? 30 : 1;
    473532
    474533/* Does this cell look kosher? */
    475534         if(OTA[cell].rms*OTA[cell].rms >
    476             OTA[cell].sky/MIN_EADU + MAX_READ_NOISE*MAX_READ_NOISE) {
     535            OTA[cell].sky/MIN_EADU + i*i*MAX_READ_NOISE*MAX_READ_NOISE) {
    477536            if(VERBOSE > 0) {
    478537               fprintf(stderr, "logonly: cell %d is unreasonably noisy, skipping...\n", cell);
     538               fprintf(stderr, "logonly: bias = %d sky = %d rms = %d satfrac = %.2f\n",
     539                       OTA[cell].bias, OTA[cell].sky, OTA[cell].rms, OTA[cell].satfrac);
    479540            }
    480541            continue;
     
    506567
    507568         } else {
    508 
    509569/* Restore the old burns */
    510570            burn_restore(naxis1-ovrscan1, naxis2-ovrscan2, naxis1,
     
    514574/* Write the corrected data back to the FITS. */
    515575         if(update) {
    516            if (ofilename) {
    517              fh_ehu(ehu, 0);    /* Seek back to the start of data */
    518 
    519              if (fh_write(ehu, outfd)) {
    520                fprintf(stderr, "\rerror: Trouble writing EXT %d to %s\n", ext, ofilename);
    521                exit(EXIT_FAILURE);
    522              }
    523 
    524              if (fh_write_padded_image(ehu, outfd, buf,
    525                                        naxis1*naxis2*naxis3*sizeof(short),
    526                                        FH_TYPESIZE_16) != FH_SUCCESS) {
     576            fh_ehu(ehu, 0);     /* Seek back to the start of data */
     577            if (fh_write_padded_image(ehu, fh_file_desc(ehu), buf,
     578                                      naxis1*naxis2*naxis3*sizeof(short),
     579                                      FH_TYPESIZE_16) != FH_SUCCESS) {
    527580               fprintf(stderr, "\rerror: failed to re-write image data for extension `%s'.\n",
    528581                       extname);
    529582               free(buf);
    530583               exit(EXIT_FAILURE);
    531              }
    532            } else {
    533              fh_ehu(ehu, 0);    /* Seek back to the start of data */
    534              if (fh_write_padded_image(ehu, fh_file_desc(ehu), buf,
    535                                        naxis1*naxis2*naxis3*sizeof(short),
    536                                        FH_TYPESIZE_16) != FH_SUCCESS) {
    537                fprintf(stderr, "\rerror: failed to re-write image data for extension `%s'.\n",
    538                        extname);
    539                free(buf);
    540                exit(EXIT_FAILURE);
    541              }
    542            }
    543          }
    544       }
    545      
     584            }
     585         }
     586      }
     587
    546588      free(buf);
    547589   }
    548590
    549    if (ofilename) {
    550      if (close(outfd)) {
    551        fprintf(stderr, "\rerror: trouble writing data to output file %s\n", ofilename);
    552        exit(EXIT_FAILURE);
    553      }
    554    }
    555      
    556591/* Dump out the postage stamp file */
    557592   if(psffile != NULL) {
     
    564599   }
    565600
     601/* Write burn info to FITS file. */
     602   if(update) persist_fits_write(OTA, ihu);
     603
     604   
     605   if(restore) {   
     606     /* Indicate in the header that the burns are not applied. */
     607     fh_set_bool(ihu, FH_AUTO, PHU_NAME_BURN_APPLIED,
     608                 FH_FALSE, PHU_COMMENT_BURN_APPLIED);
     609     fh_rewrite(ihu);
     610   }
     611
    566612   fh_destroy(ihu);
    567613
    568614/* Write the persistence data for the next image */
    569615   if(burnfile != NULL) persist_write(OTA, burnfile);
     616
    570617
    571618   exit(EXIT_SUCCESS);
     
    613660                      CELL *cell)
    614661{
    615    int i, j, k, n;
     662   int i, j, k, n, nsat;
    616663
    617664/* Get bias stats */
     
    624671
    625672/* Get sky stats */
    626    for(k=n=0; k<nx*ny; k+=((617*nx)/1000)) {
     673   for(k=n=nsat=0; k<nx*ny; k+=((617*nx)/1000)) {
    627674      i = k % nx;
    628675      j = k / nx;
    629       median_buf[n++] = data[i+NX*j];
    630    }
     676      if(data[i+NX*j] > SAT4SURE) {
     677         nsat++;
     678      } else if(data[i+NX*j] != NODATA) {
     679         median_buf[n++] = data[i+NX*j];
     680      }
     681   }
     682
     683   if(n < 20) {         /* Better have hit at least 20! */
     684      cell->sky = cell->rms = 0.0;
     685      return(-1);
     686   }
     687   cell->satfrac = ((double)nsat) / n;
     688
    631689/* First pass at sky and quartile */
    632690   cell->sky = int_median(n, median_buf);
     
    656714   printf(" restore={t|f}  Restore the input MEF by adding input fits?\n");
    657715   printf(" in=fname       Input file for previous burn persistence streaks\n");
     716   printf(" infits=fname   Input FITS file for previous burn persistence streaks (stored\n");
     717   printf("                   in table extensions).  If both this and the 'in' input\n");
     718   printf("                   file option are specified, then 'in' takes precedence.\n");
    658719   printf(" out=fname      Output file for burn streaks\n");
    659    printf(" outimage=imname Output file for MEF image\n");
     720
    660721//   printf(" trailin=fname  Input file for previous burn persistence streaks\n");
     722//   printf(" trailinfits=fname  Input FITS file for previous burn persistence streaks\n");
    661723//   printf(" trailout=fname Output file for burn streaks\n");
     724   printf(" deltables=fname Copy mef_file to new FITS file 'fname' with burn streak\n");
     725   printf("                   tables removed. NOTE: if specified, all other options\n");
     726   printf("                   will be ignored!\n");
    662727   printf(" psf=fname      Output file for PSF FITS stamp gallery\n");
    663728   printf(" psfstat=fname  Output file for PSF statistics listing\n");
     
    676741   printf(" rmask=X        Diameter growth factor of burned spots\n");
    677742   printf(" bmask=X        Box size growth of burn/star boxes\n");
     743   printf(" expire=N       Retire a blasted burn after N seconds\n");
    678744   printf(" quiet={t|f}    Quiet?\n");
    679745   printf(" verbose=N      Set verbosity bits:\n");
  • trunk/extsrc/gpcsw/gpcsrc/fits/burntool/burntool.h

    r23490 r23924  
     1#ifndef _INCLUDED_burntool_
     2#define _INCLUDED_burntool_
     3
    14/* Burn correction routines:
    25   ------------------------- 
     
    6366
    6467#define MAXCELL 64              /* Max cells in an OTA */
    65 #define MAXSIZE 700             /* Maximum vertical cell size */
     68#define MAXSIZE 2048            /* Maximum vertical cell size */
    6669
    6770#define STAR_RADIUS  4          /* Radius over which a star ctr must be max */
     
    8790#define BURN_PWR  1             /* Power law */
    8891#define BURN_EXP  2             /* Exponential */
    89 #define PSF_STAR  3             /* Unfitted: good psf star */
     92#define BURN_BLASTED 3          /* Blasted top to bottom: flag only for IPP */
     93#define PSF_STAR  9             /* Unfitted: good psf star */
    9094
    9195/* Fit error codes */
     
    9397#define FIT_TOP_ERROR  2        /* Saturation extends to top: no points */
    9498#define FIT_SLOPE_ERROR  3      /* Unreasonable fit */
     99#define FIT_EXPIRED 9           /* Don't carry any more as persistent */
    95100
    96101/* Fit parameters */
     
    152157      int rms;                  /* RMS in the sky */
    153158      int time;                 /* PON time of this cell */
     159      double satfrac;           /* Fraction of SAT4SURE saturated pixels */
    154160      int nburn;                /* Number of trails left by sat stars */
    155161      OBJBOX *burn;             /* Stars which we think left a trail */
     
    212218int write_2ddata(int nx, int ny, int *ntot, IMTYPE *data, int fd);
    213219int write_3dend(int *ntot, int fd);
     220
     221#endif /* _INCLUDED_burntool_ */
  • trunk/extsrc/gpcsw/gpcsrc/fits/burntool/man/burntool.1

    r23490 r23924  
    5353        power laws; downward burns as exponentials.
    5454       
     55        Burntool also identifies really blasted areas which are saturated from
     56        top to bottom.  These cannot be fitted, but are carried along for
     57        EXPIRE_TRAIL_TIME seconds (default 2000) in the persistence list.  The
     58        idea is to give IPP a clue that something bad happened, and of course
     59        the persistence info tells when it happened.  This expiration time
     60        for carrying the region can be set with "expire=N".
     61
    5562        Burntool finishes by rewriting the modified MEF and all of the fits in
    5663        an output file ("out=file_name").  If desired, burntool can be run
     
    112119           npsf=3          Number of contributing PSF stars
    113120           fwhm=5.38       PSF FWHM [pix] (0.0 if fit fails)
    114            m2=11.08        Second moment of distribution [pix^2]
    115            qp=5.154        Plus quadrupole qxx-qyy [pix^2]       
    116            qc=0.560        Cross quadrupole 2*qxy [pix^2]       
    117            qt=-5.201       Tangential quadrupole [pix^2]         
     121           fmed=5.76       FWHM averaged over psfavg (0.0 if fit fails)
     122           m2=11.08        Second moment of distribution [pix^2] (-99.99 if fails)
     123           qp=5.154        Plus quadrupole qxx-qyy [pix^2] (-99.99 if fails)
     124           qc=0.560        Cross quadrupole 2*qxy [pix^2] (-99.99 if fails)
     125           qt=-5.201       Tangential quadrupole [pix^2] (-99.99 if fails)
    118126                                q+ * cos(2*phi) + qx * sin(2*phi)
     127           qpm=5.154       Plus quadrupole averaged over psfavg
     128           qcm=0.560       Cross quadrupole averaged over psfavg
     129           qtm=-5.201      Tangential averaged over psfavg
    119130
    120131        Helpful utilities include:
     
    148159                Input file for previous burn persistence streaks
    149160
     161        infits=fname       
     162                Input FITS file for previous burn persistence streaks
     163                (burn streaks stored in table extensions). If both
     164                this and the 'in' input file option are specified, then
     165                'in' takes precedence.
     166
    150167        out=fname     
    151168                Output file for burn streaks
     169
     170        deltables=fname
     171                Copy mef_file to new FITS file 'fname' with burn streak
     172                tables removed. NOTE: if specified, all other options
     173                will be ignored!
    152174
    153175        psf=fname     
     
    199221        bmask=X       
    200222                Box size growth of burn/star (BMASK_GROW)
     223
     224        expire=N
     225                Retire a blasted burn after N seconds
    201226
    202227        quiet={t|f}   
  • trunk/extsrc/gpcsw/gpcsrc/fits/burntool/persistfix.c

    r23490 r23924  
    2727/* Fix up all the persistence streaks */
    2828   for(k=0; k<cell->npersist; k++) {
     29
     30/* Is this just a blasted area being carried for IPP? */
     31      if((cell->persist)[k].func == BURN_BLASTED) {
     32         if(cell->time - (cell->persist)[k].time > EXPIRE_TRAIL_TIME) {
     33            (cell->persist)[k].fiterr = FIT_EXPIRED;
     34         }
     35         continue;
     36      }
    2937
    3038/* Fit the trail */
  • trunk/extsrc/gpcsw/gpcsrc/fits/burntool/persistio.c

    r23490 r23924  
    233233   fprintf(fp,  "# Cell: %d  sky= %d   rms= %d   bias= %d\n",
    234234           cell[0].cell, cell[0].sky, cell[0].rms, cell[0].bias);
    235    fprintf(fp, "#Cell time    cx  cy  max   sx  sy   ex  ey  y0m y0p y1m y1p x0m x0p x1m x1p     slope    zero\n");
     235   fprintf(fp, "#Cell time    cx  cy  max   y0   sx  sy   ex  ey  y0m y0p y1m y1p x0m x0p x1m x1p F up    slope nfit sxf exf\n");
    236236
    237237   for(j=0; j<MAXCELL; j++) {
  • trunk/extsrc/gpcsw/gpcsrc/fits/burntool/psfstats.c

    r23525 r23924  
    1 /* psfstats.c - calculate some statistics about a PSF */
     1/* dummy version of psfstats.c - disable use of the psf */
    22#include <stdio.h>
    33#include <stdlib.h>
     
    1212#include "psf/psf.h"
    1313
    14 //#define PSFSTATS              /* ignore psf stats so we can skip fortan */
    15 //#define TEST          /* Dump all the PSF info? */
    16 
    17 static unsigned short int *psfbuf=NULL;
    18 static int npsfbuf=0;
    19 #define PHONYSKY 1000
    20 
    2114/****************************************************************/
    2215/* psf_stats(): Tell us about this PSF star */
     
    2417                     double *fwhm, double *q)
    2518{
    26 
    27 #ifdef PSFSTATS
    28 
    29    int i, j, r, r2, err, i0, j0;
    30    double d, qxx, qxy, qyy, mx, my, sum;
    31    PSF_PARAM psfout;      /* Results of fit */
    32    PSF_EXTRA psfextra;    /* Extra results from fit */
    33    double pi=4*atan(1.0);
    34 
    35    if(nx*ny > npsfbuf) {
    36       if(psfbuf != NULL) free(psfbuf);
    37       psfbuf = (ushort *)calloc(nx*ny, sizeof(ushort));
    38       npsfbuf = nx*ny;
    39    }
    40    for(i=0; i<nx*ny; i++) {
    41       j = ((int)data[i]) - bias + PHONYSKY;
    42       if(j >= 0 && j < 65535) psfbuf[i] = j;
    43       else psfbuf[i] = 0;
    44    }
    45 
    46 /* 2D PSF fit */
    47    err = psf(nx, ny, psfbuf, PSF_2DIM, NULL, &psfout, &psfextra);
    48 
    49 #ifdef TEST
    50    printf("err = %d\n", err);
    51 #endif
    52 
    53    if(!err) {
    54 #ifdef TEST
    55       printf("err = %d\n", err);
    56       printf("bias = %d\n", bias);
    57       printf("ix, iy = %d %d\n", psfout.ix, psfout.iy);
    58       printf("x0, y0 = %.2f %.2f\n", psfout.x0, psfout.y0);
    59       printf("fwhm = %.2f\n", psfout.fwhm);
    60       printf("peak, bckgnd = %.1f %.1f\n", psfout.peak, psfout.bkgnd);
    61       printf("flux, S/N = %.1f %.2f\n", psfout.flux, psfout.sn);
    62       printf("weight = %.2f\n", psfout.weight);
    63       printf("FWmaj, FWmin = %.2f %.2f\n", psfextra.majfw, psfextra.minfw);
    64       printf("theta = %.1f\n", psfextra.thfw*180/pi);
    65       printf("wpeak, wbkgnd = %.1f %.1f\n", psfextra.wpeak, psfextra.wbkgnd);
    66       printf("dflux, dbckgnd = %.1f %.1f\n", psfextra.dflux, psfextra.dbkgnd);
    67       printf("rmsbkgnd = %.1f\n", psfextra.rmsbkgnd);
    68 #endif
    69       fwhm[0] = psfextra.majfw;
    70       fwhm[1] = psfextra.minfw;
    71       fwhm[2] = psfextra.thfw;
    72 
    73 /* Moments */
    74       i0 = MIN(psfout.ix, nx-psfout.ix);
    75       j0 = MIN(psfout.iy, ny-psfout.iy);
    76    } else {
    77       fwhm[0] = fwhm[1] = fwhm[2] = 0.0;
    78       i0 = nx/2;
    79       j0 = ny/2;
    80    }
    81 
    82    r = MIN(i0, j0) - 1;
    83    sum = mx = my = qxx = qxy = qyy = 0.0;
    84    for(j=j0-r; j<j0+r; j++) {
    85       for(i=i0-r; i<i0+r; i++) {
    86          r2 = (i-i0)*(i-i0) + (j-j0)*(j-j0);
    87          if(r2 > r*r) continue;
    88          d = (double)data[i+j*nx] - bias;
    89          sum += d;
    90          mx += (i+0.5) * d;
    91          my += (j+0.5) * d;
    92          qxx += (i+0.5) * (i+0.5) * d;
    93          qxy += (i+0.5) * (j+0.5) * d;
    94          qyy += (j+0.5) * (j+0.5) * d;
    95       }
    96    }
    97    if(sum > 0) {
    98       mx /= sum;
    99       my /= sum;
    100       qxx = qxx / sum - mx*mx;
    101       qxy = qxy / sum - mx*my;
    102       qyy = qyy / sum - my*my;
    103    } else {
    104       q[0] = q[1] = q[2] = 0.0;
    105    }
    106 
    107    q[0] = (qxx + qyy) / 2;
    108    q[1] = qxx - qyy;
    109    q[2] = 2 * qxy;
    110 
    111 /* phi is the angular position of the OTA around the center [rad] */
    112 //   qtang = (qxx - qyy) * cos(2*phi) + 2 * qxy * sin(2*phi);
    113 
    114 #ifdef TEST
    115    printf("Sum= %.1f mx= %.2f my= %.2f qxx= %.4f qxy= %.4f qyy= %.4f\n",
    116           sum, mx, my, qxx, qxy, qyy);
    117 #endif
    118 
    119 # else
    120    printf("psfstats ignored in this build\n");
    121    return (1);
    122 #endif /* PSFSTATS */
    123    return(0);
     19   return(1);
    12420}
  • trunk/extsrc/gpcsw/gpcsrc/fits/burntool/trailfit.c

    r23490 r23924  
    148148   }
    149149
     150/* Burn extends all the way to the top */
    150151   if(up && y1 >= ny-1) {
    151152      box->slope = -1.0;
    152       box->nfit = xe - xs + 1;
    153153      box->fiterr = FIT_TOP_ERROR;
     154/* No fit here because the burn is to top, but the persist may be fitable */
     155      if(box->sy > FIT_EDGE) {
     156         box->nfit = xe - xs + 1;
     157/* However, if burn extends to bottom, persist will not be able to fit */
     158      } else {
     159         box->nfit = 1;         /* Just enough to keep it alive */
     160         box->func = BURN_BLASTED;
     161      }
    154162      return(0);
    155163   }
  • trunk/extsrc/gpcsw/gpcsrc/fits/libfh/Makefile

    r23490 r23924  
    1010include ../Make.Common
    1111
    12 VERSION = 2.01
     12VERSION = 2.03
    1313CCWARN += $(WERROR)
    1414CCDEFS += -DHAVE_FH_VALIDATE
    15 SRCS = fh.c
     15SRCS = fh.c fh_table.c
    1616HDRS = fh.h fh_registry.h fh_registry.asm
    1717
  • trunk/extsrc/gpcsw/gpcsrc/fits/libfh/fh.h

    r23490 r23924  
    3636#define FH_COUNT_CARDS "@FH_COUNT_CARDS@"  /* See also: fh_count_cards() */
    3737
     38/* Size of buffer in fhTableCol structure for format string. */
     39#define FH_TABLE_FORMAT_STR_LEN 20
     40
    3841typedef enum /* fh_result -- Result Code for most functions in this library: */
    3942{  FH_SUCCESS = 0,
     
    6467} fh_bool;                      /* Value for "logical" FITS cards */
    6568
     69/* The possible field formats in an ASCII table extension. */
     70typedef enum
     71{
     72   FH_TABLE_FORMAT_CHAR = 'A',   /* Designates a string field. */
     73   FH_TABLE_FORMAT_INT = 'I',    /* Integer field. */
     74   FH_TABLE_FORMAT_FLOAT = 'F',  /* Float field. */
     75   FH_TABLE_FORMAT_DOUBLE = 'D'  /* Double field. */
     76} fhTableFormat;
     77
     78typedef struct
     79{
     80  char * name;    /* Short string identifying this column's data. */
     81  char * comment; /* Longer blurb about this column's data. */
     82  char * units;   /* Units for column's data. */
     83  fhTableFormat format; /* Type of data held in this column. */
     84  int width;      /* Width of this column's data, in characters. */
     85  int dec_places; /* For doubles, etc, the number of decimal places to write. */
     86  /* Everything below this line is generated by fh_table_init(). These
     87   * do not need to be filled out by hand. */
     88  int first_char_pos; /* Position of first character in this column. Calculated on the fly. */
     89  char format_string[FH_TABLE_FORMAT_STR_LEN]; /* Format string for printf. */
     90} fhTableCol;
     91
     92/* Describes the layout of an ASCII table. */
     93typedef struct
     94{
     95  char * extname;    /* Name of table extension. */
     96  int num_cols;      /* Number of rows in table. */
     97  int num_rows;      /* Number of columns in table. */
     98  fhTableCol * cols; /* Array of cols columns. */
     99  /* Everything below this line is generated by fh_table_init(). These
     100   * do not need to be filled out by hand. */
     101  int widest_col;    /* Width of widest column. */
     102  char * strbuf;     /* Will point to buf of widest_col + 1 chars. */
     103  int row_width;     /* Number of characters in a table tow. */
     104  int table_size;    /* Total size of table data, in bytes. */
     105} fhTable;
     106
    66107typedef void* HeaderUnit; /* Handle to a list of FITS cards allocated by fh_create. */
    67108
     
    345386double fh_idx(HeaderUnit hu); /* idx of the last card returned by fh_next */
    346387fh_result fh_merge(HeaderUnit hu, const HeaderUnit source); /* source unchanged */
     388
     389/* ---------------------------------------------------------
     390 * Reading and writing tables
     391 * ---------------------------------------------------------
     392 */
     393
     394fh_result
     395fh_table_init(fhTable * table);
     396/* Call once to initialise information about a given table. This will
     397 * walk through all of the columns and fill out starting character
     398 * positions and so forth, and will calculate the size of each row
     399 * and the size of the table as a whole, and so forth. */
     400
     401fh_result
     402fh_table_populate_header(HeaderUnit hu, fhTable * table);
     403/* Sets up keywords for a table extension header, including all
     404 * of the column information, etc. This will initialise the table
     405 * structure if necessary. This does not include any fh_write...
     406 * calls to write the header to the file. */
     407
     408fh_result
     409fh_reserve_padded_table(HeaderUnit hu, int fd);
     410/* Analogous to fh_reserve_padded_image() for tables. This reserves
     411 * space for the table based on the information in the passed extension
     412 * header at the current location in the file. The header must therefore
     413 * be populated prior to calling this function (use
     414 * fh_table_populate_header() to build a header from a table structure.
     415 */
     416
     417fh_result
     418fh_map_table(HeaderUnit hu, void** data, int size);
     419fh_result
     420fh_munmap_table(HeaderUnit hu);
     421/* Analogous to fh_map_raw_image() and fh_munmap_image() for tables.
     422 * These calls give the calling program access to a memory mapped pointer
     423 * to the data in the FITS file.  The data format is equivalent to what
     424 * would be produced by fh_read_image() with typesize=FH_TYPESIZE_RAW.
     425 * Whether or not the file can be modified through the memory-mapped
     426 * pointer depends on whether fh_mode=FH_FILE_RDONLY or FH_FILE_RDWR
     427 * in the call to fh_file() that opened the file.  Be sure to use
     428 * fh_munmap_table() when done.
     429 *
     430 * `size' must match the exact number bytes of data which should be in
     431 * the file, according to its NAXIS and BITPIX values, excluding padding.
     432 *
     433 * %%% TODO: FH_FILE_RDWR is completely unimplemented!
     434 *
     435 * %%% TODO: File locking is not handled properly yet.  For now, use
     436 * FH_FILE_RDONLY_NOLOCK and FH_FILE_RDWR_NOLOCK, especially on MEF
     437 * files where the locks are left in place until the file is closed
     438 * otherwise!
     439 */
     440
     441fh_result
     442fh_table_read_value(fhTable * table, void * data,
     443                    int row, int col, void * value);
     444/* Reads a single value from an existing table at (row, column).
     445 * The table field is interpreted based on the format of that
     446 * column in the table structure, and is passed back via the
     447 * 'value' pointer. It is the caller's responsibility to ensure
     448 * that 'value' points to a location with sufficient space for
     449 * the data in the field, and that any casting to & from a void
     450 * pointer is done correctly.
     451 */
     452
     453fh_result
     454fh_table_write_value(fhTable * table, void * data,
     455                     int row, int col, void * value);
     456/* Similar to read_value() above, this function overwrites a field
     457 * in a table at the given location.  The caller is repsonsible
     458 * for casting the location of the stored data to a void pointer.
     459 */
     460
     461
    347462#endif /* _INCLUDED_fh */
  • trunk/extsrc/gpcsw/gpcsrc/fits/libfhreg/gpc_telescope.h

    r23490 r23924  
    8484MK_PFL(  522.5, ROT     ,6, "Telescope rotator angle (degrees)"               )
    8585MK_PFL(  522.6, POSANGLE,6, "Telescope position angle (degrees)"              )
    86 MK_PFL(  522.7, PA_DEBUG,6, "Telescope position angle (alternate value)"      )
    8786MK_PFL(  523.0, COMRA   ,6, "Commanded telescope Right Ascension (degrees)"   )
    8887MK_PFL(  523.1, COMDEC  ,6, "Commanded telescope Declination (degrees)"       )
Note: See TracChangeset for help on using the changeset viewer.