IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
May 22, 2013, 4:54:28 PM (13 years ago)
Author:
eugene
Message:

fiximids seems to be working correctly now

Location:
branches/eam_branches/ipp-20130509/Ohana/src/uniphot
Files:
9 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20130509/Ohana/src/uniphot/Makefile

    r35103 r35593  
    1 default: uniphot setphot setphot_client setfwhm
     1default: uniphot setphot setphot_client setfwhm setposangle setposangle_client setastrom setastrom_client fiximids fiximids_client
    22help:
    33@echo "make options: uniphot setphot setfwhm default help install default (uniphot setphot setfwhm)"
    … …  
    2525setastrom: $(BIN)/setastrom.$(ARCH)
    2626setastrom_client: $(BIN)/setastrom_client.$(ARCH)
     27fiximids: $(BIN)/fiximids.$(ARCH)
     28fiximids_client: $(BIN)/fiximids_client.$(ARCH)
    2729
    28 install: $(DESTBIN)/uniphot $(DESTBIN)/setfwhm $(DESTBIN)/setphot $(DESTBIN)/setphot_client $(DESTBIN)/setposangle $(DESTBIN)/setposangle_client $(DESTBIN)/setastrom $(DESTBIN)/setastrom_client
     30install: $(DESTBIN)/uniphot $(DESTBIN)/setfwhm $(DESTBIN)/setphot $(DESTBIN)/setphot_client $(DESTBIN)/setposangle $(DESTBIN)/setposangle_client $(DESTBIN)/setastrom $(DESTBIN)/setastrom_client $(DESTBIN)/fiximids $(DESTBIN)/fiximids_client
    2931
    3032UNIPHOT = \
    … …  
    137139$(SETFWHM): $(INC)/uniphot.h
    138140$(BIN)/setfwhm.$(ARCH): $(SETFWHM)
     141
     142FIXIMIDS = \
     143$(SRC)/fiximids.$(ARCH).o                \
     144$(SRC)/initialize_fiximids.$(ARCH).o     \
     145$(SRC)/load_images_fiximids.$(ARCH).o    \
     146$(SRC)/ImageSubsetFixImageIDs.$(ARCH).o  \
     147$(SRC)/update_dvo_fiximids.$(ARCH).o     \
     148$(SRC)/update_catalog_fiximids.$(ARCH).o \
     149$(SRC)/SetSignals.$(ARCH).o              \
     150$(SRC)/Shutdown.$(ARCH).o           
     151
     152$(FIXIMIDS): $(INC)/fiximids.h
     153$(BIN)/fiximids.$(ARCH): $(FIXIMIDS)
     154
     155FIXIMIDS_CLIENT = \
     156$(SRC)/fiximids_client.$(ARCH).o         \
     157$(SRC)/initialize_fiximids.$(ARCH).o     \
     158$(SRC)/load_images_fiximids.$(ARCH).o    \
     159$(SRC)/ImageSubsetFixImageIDs.$(ARCH).o  \
     160$(SRC)/update_dvo_fiximids.$(ARCH).o     \
     161$(SRC)/update_catalog_fiximids.$(ARCH).o \
     162$(SRC)/SetSignals.$(ARCH).o              \
     163$(SRC)/Shutdown.$(ARCH).o
     164
     165$(FIXIMIDS_CLIENT): $(INC)/fiximids.h
     166$(BIN)/fiximids_client.$(ARCH): $(FIXIMIDS_CLIENT)
     167

  • branches/eam_branches/ipp-20130509/Ohana/src/uniphot/include/fiximids.h

    r35591 r35593  
    22# include <dvo.h>
    33# include <signal.h>
    4 
    5 typedef struct {
    6   unsigned int start;
    7   unsigned int stop;
    8   Coords coords;
    9 } Mosaic;
    104
    115// to determine the image ID, we need time (tzero + trate) and photcode
    … …  
    148  unsigned short   NX;                   // image width
    159  unsigned short   NY;                   // image height
     10  unsigned int     nstar;                // number of stars on the image
    1611  e_time           tzero;                // readout time (row 0)
    1712  unsigned char    trate;                // scan rate (100 usec/pixel)
    … …  
    2823char        *SINGLE_CPT;
    2924int          VERBOSE;
    30 int          RESET;
     25int          VERBOSE_IMSTATS;
     26int          SUMMARY_IMSTATS;
     27char        *IMSTATS_FILE; // location of partial imstats file for remote client
    3128int          UPDATE;
    3229int          PARALLEL;
    … …  
    4037void          ConfigInit                        PROTO((int *argc, char **argv));
    4138
    42 void          usage_fiximids                 PROTO(());
    43 void          initialize_fiximids            PROTO((int argc, char **argv));
    44 int           args_fiximids                  PROTO((int argc, char **argv));
    45 
    46 void          usage_fiximids_client          PROTO(());
    47 void          initialize_fiximids_client     PROTO((int argc, char **argv));
    48 int           args_fiximids_client           PROTO((int argc, char **argv));
    49 
    50 // int           update_dvo_fiximids_client  PROTO((ImageSubset *image, off_t Nimage, FlatCorrectionTable *flatcorr));
     39void          usage_fiximids                    PROTO(());
     40void          initialize_fiximids               PROTO((int argc, char **argv));
     41int           args_fiximids                     PROTO((int argc, char **argv));
     42                                               
     43void          usage_fiximids_client             PROTO(());
     44void          initialize_fiximids_client        PROTO((int argc, char **argv));
     45int           args_fiximids_client              PROTO((int argc, char **argv));
    5146
    5247ImageSubset  *ImageSubsetLoad                   PROTO((char *filename, off_t *nimage));
    5348int           ImageSubsetSave                   PROTO((char *filename, ImageSubset *image, off_t Nimage));
    54 Image        *ImagesFromSubset                  PROTO((ImageSubset *subset, off_t N));
    5549ImageSubset  *ImagesToSubset                    PROTO((Image *image, off_t N));
    5650
    57 void          lock_image_db                     PROTO((FITS_DB *db, char *filename));
    58 void          unlock_image_db                   PROTO((FITS_DB *db));
    59 void          create_image_db                   PROTO((FITS_DB *db));
    6051void          set_db                            PROTO((FITS_DB *in));
    6152int           Shutdown                          PROTO((char *format, ...) OHANA_FORMAT(printf, 1, 2) );
    … …  
    6455int           SetSignals                        PROTO((void));
    6556
    66 Image        *load_images_fiximids           PROTO((FITS_DB *db, off_t *Nimage));
    67 int           update_dvo_fiximids            PROTO((ImageSubset *image, off_t Nimage));
    68 int           update_dvo_fiximids_parallel   PROTO((SkyTable *sky, ImageSubset *image, off_t Nimage));
     57Image        *load_images_fiximids              PROTO((FITS_DB *db, off_t *Nimage));
     58int           update_dvo_fiximids               PROTO((ImageSubset *image, off_t Nimage));
     59int           update_dvo_fiximids_parallel      PROTO((SkyTable *sky, ImageSubset *image, off_t Nimage));
     60                                               
     61void          update_catalog_fiximids           PROTO((Catalog *catalog));
    6962
    70 void          update_catalog_fiximids        PROTO((Catalog *catalog, ImageSubset *image, off_t *index, off_t Nimage));
    71 int           fiximids_local_astrometry      PROTO((float *posAngle, float *pltScale, double x, double y, Coords *mosaic, Coords *coords));
     63void          sort_image_times                  PROTO((e_time *S, e_time *E, short *C, unsigned int *I, unsigned int *Q, off_t N));
     64int           FindImageID                       PROTO((off_t *ID, off_t *Seq, e_time time, short photcode));
     65void          initImageIndex                    PROTO((ImageSubset *image, off_t Nimage_init));
    7266
    73 off_t         getMosaicByTimes                  PROTO((unsigned int start, unsigned int stop, unsigned int *startMos, unsigned int *stopMos, off_t *indexMos));
    74 void          sort_mosaic_times                 PROTO((unsigned int *S, unsigned int *E, off_t *I, off_t N));
    75 void          initMosaics                       PROTO((ImageSubset *image, off_t Nimage));
    76 Mosaic       *getMosaicForImage                 PROTO((off_t im));
     67void          BumpValidImage                    PROTO((int Seq));
     68void          BumpInvalidImage                  PROTO((int Seq));
     69void          CompareImageCounts                PROTO((ImageSubset *image, off_t Nimage_comp));
     70void          SummaryImageStats                 PROTO((ImageSubset *image, off_t Nimage_comp));
     71
     72int           ImageValidSave                    PROTO((char *filename));
     73int           ImageValidLoad                    PROTO((char *filename));

  • branches/eam_branches/ipp-20130509/Ohana/src/uniphot/src/ImageSubsetFixImageIDs.c

    r35591 r35593  
    11# include "fiximids.h"
     2# define FT_BZERO_INT16 1.0*0x8000       
     3# define FT_BZERO_INT32 1.0*0x80000000
    24
    35# define GET_COLUMN(OUT,NAME,TYPE) \
    … …  
    810
    911  int Ncol;
    10   off_t i, j, k;
     12  off_t i;
    1113  off_t Nrow;
    1214  Header header;
    … …  
    5557
    5658  GET_COLUMN (imageID,    "IMAGE_ID",   int);
    57   GET_COLUMN (NX,         "NX",         int);
    58   GET_COLUMN (NY,         "NY",         int);
    59   GET_COLUMN (tzero,      "TZERO",      float);
    60   GET_COLUMN (trate,      "TRATE",      float);
     59  GET_COLUMN (NX,         "NX",         short);
     60  GET_COLUMN (NY,         "NY",         short);
     61  GET_COLUMN (nstar,      "nstar",      int);
     62  GET_COLUMN (tzero,      "TZERO",      int);
     63  GET_COLUMN (trate,      "TRATE",      short);
    6164  GET_COLUMN (photcode,   "PHOTCODE"  , short);
    6265
    … …  
    6669    image[i].NX                = NX[i];
    6770    image[i].NY                = NY[i];
     71    image[i].nstar             = nstar[i];
    6872    image[i].tzero             = tzero[i];
    6973    image[i].trate             = trate[i];
    … …  
    7579  free (NX);
    7680  free (NY);
     81  free (nstar);
    7782  free (tzero);
    7883  free (trate);
    … …  
    8590int ImageSubsetSave(char *filename, ImageSubset *image, off_t Nimage) {
    8691
    87   off_t i, j, k;
     92  off_t i;
    8893  Header header;
    8994  Header theader;
    … …  
    100105  // an unsigned int needs to have bzero of 0x8000
    101106  gfits_define_bintable_column (&theader, "J",   "IMAGE_ID",   "image ID", "tmp", 1.0, FT_BZERO_INT32);
    102   gfits_define_bintable_column (&theader, "J",   "NX",         "tmp",      "tmp", 1.0, 0.0);
    103   gfits_define_bintable_column (&theader, "J",   "NY",         "tmp",      "tmp", 1.0, 0.0);
    104   gfits_define_bintable_column (&theader, "E",   "TZERO",      "tmp",      "tmp", 1.0, 0.0);
    105   gfits_define_bintable_column (&theader, "E",   "TRATE",      "tmp",      "tmp", 1.0, 0.0);
     107  gfits_define_bintable_column (&theader, "I",   "NX",         "tmp",      "tmp", 1.0, 0.0);
     108  gfits_define_bintable_column (&theader, "I",   "NY",         "tmp",      "tmp", 1.0, 0.0);
     109  gfits_define_bintable_column (&theader, "J",   "NSTAR",      "tmp",      "tmp", 1.0, FT_BZERO_INT32);
     110  gfits_define_bintable_column (&theader, "J",   "TZERO",      "tmp",      "tmp", 1.0, 0.0);
     111  gfits_define_bintable_column (&theader, "I",   "TRATE",      "tmp",      "tmp", 1.0, 0.0);
    106112  gfits_define_bintable_column (&theader, "I",   "PHOTCODE",   "tmp",      "tmp", 1.0, 0.0);
    107113
    … …  
    109115  gfits_create_table (&theader, &ftable);
    110116
    111   unsigned int *imageID;
    112   int *NX, *NY;
    113   float *tzero, *trate;
     117  unsigned int *imageID, *nstar;
     118  short *NX, *NY;
     119  e_time *tzero;
     120  short *trate;
    114121  short *photcode;
    115122
    116123  // create intermediate storage arrays
    117124  ALLOCATE (imageID,    unsigned int, Nimage);
    118   ALLOCATE (NX,         int,          Nimage);
    119   ALLOCATE (NY,         int,          Nimage);
    120   ALLOCATE (tzero,      float,        Nimage);
    121   ALLOCATE (trate,      float,        Nimage);
    122   ALLOCATE (photcode,   int,          Nimage);
     125  ALLOCATE (NX,         short,        Nimage);
     126  ALLOCATE (NY,         short,        Nimage);
     127  ALLOCATE (nstar,      unsigned int, Nimage);
     128  ALLOCATE (tzero,      e_time,       Nimage);
     129  ALLOCATE (trate,      short,        Nimage);
     130  ALLOCATE (photcode,   short,        Nimage);
    123131
    124132  // assign the storage arrays
    … …  
    127135    NX[i]       = image[i].NX      ;       
    128136    NY[i]       = image[i].NY      ;       
     137    nstar[i]    = image[i].nstar   ;       
    129138    tzero[i]    = image[i].tzero   ; 
    130139    trate[i]    = image[i].trate   ; 
    … …  
    136145  gfits_set_bintable_column (&theader, &ftable, "NX",         NX,       Nimage);
    137146  gfits_set_bintable_column (&theader, &ftable, "NY",         NY,       Nimage);
     147  gfits_set_bintable_column (&theader, &ftable, "NSTAR",      nstar,    Nimage);
    138148  gfits_set_bintable_column (&theader, &ftable, "TZERO",      tzero,    Nimage);
    139149  gfits_set_bintable_column (&theader, &ftable, "TRATE",      trate,    Nimage);
    … …  
    143153  free (NX      );
    144154  free (NY      );
     155  free (nstar   );
    145156  free (tzero   );
    146157  free (trate   );
    … …  
    169180}
    170181
    171 // XXX drop these two?
    172 
    173 # if (0)
    174 Image *ImagesFromSubset (ImageSubset *subset, off_t N) {
    175 
    176   off_t i;
    177 
    178   // we have been given an ImageSubset array, containing a reduced set of image fields
    179   // create full a Image array and save the needed values
    180   Image *image = NULL;
    181   ALLOCATE (image, Image, N);
    182   for (i = 0; i < N; i++) {
    183     image[i].imageID       = subset[i].imageID;
    184     image[i].NX            = subset[i].NX     ;
    185     image[i].NY            = subset[i].NY     ;
    186     image[i].tzero         = subset[i].tzero  ;
    187     image[i].trate         = subset[i].trate  ;
    188     image[i].coords        = subset[i].coords ;
    189   }
    190   return image;
    191 }
    192 
    193182ImageSubset *ImagesToSubset (Image *image, off_t N) {
    194183
    … …  
    201190
    202191  for (i = 0; i < N; i++) {
    203     subset[i].imageID = image[i].imageID;
    204     subset[i].NX      = image[i].NX     ;
    205     subset[i].NY      = image[i].NY     ;
    206     subset[i].tzero   = image[i].tzero  ;
    207     subset[i].trate   = image[i].trate  ;
    208     subset[i].coords  = image[i].coords ;
     192    subset[i].imageID  = image[i].imageID  ;
     193    subset[i].NX       = image[i].NX       ;
     194    subset[i].NY       = image[i].NY       ;
     195    subset[i].nstar    = image[i].nstar    ;
     196    subset[i].tzero    = image[i].tzero    ;
     197    subset[i].trate    = image[i].trate    ;
     198    subset[i].photcode = image[i].photcode ;
    209199  }
    210200  return subset;
    211201}
    212202
    213 # endif (0)
     203

  • branches/eam_branches/ipp-20130509/Ohana/src/uniphot/src/fiximids.c

    r35591 r35593  
    2424  status = update_dvo_fiximids (subset, Nimage);
    2525
     26  if (VERBOSE_IMSTATS) {
     27    CompareImageCounts (subset, Nimage);
     28  }
     29
     30  if (SUMMARY_IMSTATS) {
     31    SummaryImageStats (subset, Nimage);
     32  }
     33
    2634  if (!status) exit (1);
    2735  exit (0);

  • branches/eam_branches/ipp-20130509/Ohana/src/uniphot/src/fiximids_client.c

    r35591 r35593  
    1 # include "setposangle.h"
     1# include "fiximids.h"
    22
    3 // setposangle_client is run on a remote host and is responsible for updating the catalogs
     3// fiximids_client is run on a remote host and is responsible for updating the catalogs
    44// owned by that host. 
    55
    6 // setposangle_client -hostname (hostname) -catdir (catdir) -hostdir (hostdir) -images (images)
     6// fiximids_client -hostname (hostname) -catdir (catdir) -hostdir (hostdir) -images (images)
    77// load the SkyTable and HostTable from (catdir) [contains the host responsibilities]
    88// load the image data from (images) [contains the zp infomation to apply]
    … …  
    1313
    1414  // get configuration info, args, lockfile (set CATDIR, HOST_ID, HOSTDIR, IMAGES)
    15   initialize_setposangle_client (argc, argv);
     15  initialize_fiximids_client (argc, argv);
    1616
    1717  // load the image subset table from the specified location
    … …  
    1919  ImageSubset *subset = ImageSubsetLoad (IMAGES, &Nimage);
    2020
    21   update_dvo_setposangle (subset, Nimage);
     21  update_dvo_fiximids (subset, Nimage);
    2222
    2323  exit (0);

  • branches/eam_branches/ipp-20130509/Ohana/src/uniphot/src/initialize_fiximids.c

    r35592 r35593  
    5454  if ((N = get_argument (argc, argv, "-v"))) {
    5555    VERBOSE = TRUE;
     56    remove_argument (N, &argc, argv);
     57  }
     58
     59  SUMMARY_IMSTATS = FALSE;
     60  if ((N = get_argument (argc, argv, "-imstats"))) {
     61    SUMMARY_IMSTATS = TRUE;
     62    remove_argument (N, &argc, argv);
     63  }
     64  VERBOSE_IMSTATS = FALSE;
     65  if ((N = get_argument (argc, argv, "-vimstats"))) {
     66    SUMMARY_IMSTATS = TRUE;
     67    VERBOSE_IMSTATS = TRUE;
    5668    remove_argument (N, &argc, argv);
    5769  }
    … …  
    168180  }
    169181
     182  VERBOSE_IMSTATS = FALSE;
     183  IMSTATS_FILE = NULL;
     184  if ((N = get_argument (argc, argv, "-imstats"))) {
     185    remove_argument (N, &argc, argv);
     186    IMSTATS_FILE = strcreate (argv[N]);
     187    remove_argument (N, &argc, argv);
     188  }
     189
    170190  UPDATE = FALSE;
    171191  if ((N = get_argument (argc, argv, "-update"))) {

  • branches/eam_branches/ipp-20130509/Ohana/src/uniphot/src/load_images_fiximids.c

    r35592 r35593  
    11# include "fiximids.h"
    2 
    3 off_t Nimage;
    4 e_time *startImage;
    5 e_time *stopImage;
    6 short *photcodeImage;
    7 
    8 Image *load_images_fiximids (FITS_DB *db, off_t *N) {
     2# define FT_BZERO_INT32 1.0*0x80000000
     3
     4off_t   Nimage = 0;
     5e_time *startImage = NULL;
     6e_time *stopImage = NULL;
     7short  *photcodeImage = NULL;
     8unsigned int *imageID = NULL;
     9unsigned int *imageSeq = NULL;
     10int *Nimage_valid = NULL;
     11int *Nimage_invalid = NULL;
     12
     13Image *load_images_fiximids (FITS_DB *db, off_t *Nimage_load) {
    914
    1015  Image *image;
    … …  
    1621
    1722  /* use a vtable to keep the images to be calibrated */
    18   image = gfits_table_get_Image (&db[0].ftable, Nimage, &db[0].swapped);
     23  image = gfits_table_get_Image (&db[0].ftable, Nimage_load, &db[0].swapped);
    1924  if (!image) {
    2025      fprintf (stderr, "ERROR: failed to read images\n");
    … …  
    2227  }
    2328
    24   // XXX init lookup tables
    25 
    26   fprintf (stderr, "loaded "OFF_T_FMT" images\n", *Nimage);
     29  fprintf (stderr, "loaded "OFF_T_FMT" images\n", *Nimage_load);
    2730
    2831  return (image);
     32}
     33
     34// find mosaic frames (unique time periods & photcode name matches mosaic)
     35void initImageIndex (ImageSubset *image, off_t Nimage_init) {
     36
     37  off_t i;
     38
     39  Nimage = Nimage_init;
     40  ALLOCATE (startImage, e_time, Nimage);
     41  ALLOCATE (stopImage, e_time, Nimage);
     42  ALLOCATE (photcodeImage, short, Nimage);
     43  ALLOCATE (imageID, unsigned int, Nimage);
     44  ALLOCATE (imageSeq, unsigned int, Nimage);
     45
     46  ALLOCATE (Nimage_valid, int, Nimage);
     47  ALLOCATE (Nimage_invalid, int, Nimage);
     48
     49  memset (Nimage_valid, 0, Nimage*sizeof(int));
     50  memset (Nimage_invalid, 0, Nimage*sizeof(int));
     51
     52  // save the image data in the local static arrays
     53  for (i = 0; i < Nimage; i++) {
     54
     55    /* set image time range (small boundary buffer) */
     56    e_time start = image[i].tzero - MAX(0.01*image[i].trate*image[i].NY, 1);
     57    e_time stop  = image[i].tzero + MAX(1.01*image[i].trate*image[i].NY, 1);
     58
     59    startImage[i] = start;
     60    stopImage[i] = stop;
     61    photcodeImage[i] = image[i].photcode;
     62    imageID[i] = image[i].imageID;
     63    imageSeq[i] = i;
     64  }
     65
     66  // sort the index, start, and stop by the start times:
     67  sort_image_times (startImage, stopImage, photcodeImage, imageID, imageSeq, Nimage);
     68 
     69  return;
    2970}
    3071
    … …  
    3273// This code is copied from relastro/src/MosaicOps.c
    3374
    34 off_t FindImageID (e_time time, short photcode) {
    35 
    36   // use bisection to find the overlapping imagea
    37 
    38   off_t Nlo, Nhi, N;
    39 
    40   // find the last mosaic before start
    41   Nlo = 0; Nhi = Nmosaic;
     75// find the image by time. this search makes the following assumptions:
     76// 1) there is only a single EXPOSURE from a single camera with the given exposure time
     77// 2) there may be multiple IMAGES from that exposure, distinguished by photcode
     78// 3) the image data has been sorted by start time
     79// 4) once we have (time > start[i]) we are past the images of interest
     80int FindImageID (off_t *ID, off_t *Seq, e_time time, short photcode) {
     81
     82  off_t Nlo, Nhi, N, Ni;
     83
     84  *ID = 0;
     85  *Seq = -1;
     86
     87  // Find the last image ending BEFORE 'time' (time > start[i]).
     88  // Use bisection to find the overlapping images.
     89  Nlo = 0; Nhi = Nimage;
    4290  while (Nhi - Nlo > 10) {
    4391    N = 0.5*(Nlo + Nhi);
    44     if (startMos[N] < start) {
     92    if (stopImage[N] < time) {
    4593      Nlo = MAX(N, 0);
    4694    } else {
    47       Nhi = MIN(N + 1, Nmosaic);
     95      Nhi = MIN(N + 1, Nimage);
    4896    }
    4997  }
    5098
    51   // check for the matched mosaic starting from Nlo
    52   // we may have to go much beyond Nlo since stop is not sorted
    53   // can we use a sorted version of stop to check when we are beyond the valid range??
    54   for (N = Nlo; N < Nmosaic; N++) {
    55     if (stop  < stopMos[N]) continue;
    56     if (start > startMos[N])  continue;
    57     if (stop  < startMos[N]) return (-1);
    58     return (indexMos[N]);
    59   }
    60 
    61   return (-1);
     99  // check for the matched images starting from Nlo
     100  for (N = Nlo; N < Nimage; N++) {
     101    if (time > stopImage[N]) continue;
     102    if (time < startImage[N]) return FALSE;
     103    // a possible image match (start <= time <= stop)
     104    for (Ni = N; (Ni < Nimage) && (time <= stopImage[Ni]); Ni++) {
     105      if (photcodeImage[Ni] == photcode) {
     106        *ID = imageID[Ni];
     107        *Seq = imageSeq[Ni];
     108        return TRUE;
     109      }
     110    }
     111    return FALSE;
     112  }
     113  return FALSE;
     114}
     115
     116// increment the count of valid detections for this image
     117void BumpValidImage (int Seq) {
     118  myAssert (Seq >= 0, "impossible seq");
     119  myAssert (Seq < Nimage, "impossible seq");
     120  Nimage_valid[Seq] ++;
     121}
     122
     123// increment the count of invalid detections for this image
     124void BumpInvalidImage (int Seq) {
     125  myAssert (Seq >= 0, "impossible seq");
     126  myAssert (Seq < Nimage, "impossible seq");
     127  Nimage_invalid[Seq] ++;
     128}
     129
     130void CompareImageCounts (ImageSubset *image, off_t Nimage_comp) {
     131
     132  off_t i;
     133
     134  for (i = 0; i < Nimage_comp; i++) {
     135    // XXX can I limit the images to those in range of the user region?
     136
     137    if (Nimage_valid[i] + Nimage_invalid[i] == 0) continue;     // XXX for now, skip (optionally?) images with no matched detections
     138    if (image[i].photcode == 0) continue; // skip images with 0 photcode (eg, PHU)
     139
     140    char *date = ohana_sec_to_date (image[i].tzero);
     141    char *code = GetPhotcodeNamebyCode (image[i].photcode);
     142    fprintf (stderr, "image %s %s : %d  %d  %d  : %f\n", date, code, image[i].nstar, Nimage_valid[i], Nimage_invalid[i],
     143             (Nimage_valid[i] + Nimage_invalid[i]) / (float) image[i].nstar);
     144  }
     145}
     146
     147void SummaryImageStats (ImageSubset *image, off_t Nimage_comp) {
     148
     149  off_t i;
     150
     151  // let's report the following:
     152  // Nimages with corrected image IDs
     153  // Nimages with Nsum / Nstar < 0.95
     154  // Nimages with Nsum / Nstar < 0.99
     155
     156  int Nfixed = 0;
     157  int Nbad_1 = 0;
     158  int Nbad_2 = 0;
     159
     160  for (i = 0; i < Nimage_comp; i++) {
     161    if (Nimage_valid[i] + Nimage_invalid[i] == 0) continue;     // XXX for now, skip (optionally?) images with no matched detections
     162    if (image[i].photcode == 0) continue; // skip images with 0 photcode (eg, PHU)
     163
     164    if (Nimage_invalid[i] > 0) Nfixed ++;
     165   
     166    float found_ratio = (Nimage_valid[i] + Nimage_invalid[i]) / (float) image[i].nstar;
     167   
     168    if (found_ratio < 0.99) Nbad_1 ++;
     169    if (found_ratio < 0.95) Nbad_2 ++;
     170  }
     171
     172  fprintf (stderr, "Nfixed: %d, Ndrop (>1%%) %d, Ndrop (>5%%) %d\n", Nfixed, Nbad_1, Nbad_2);
    62173}
    63174
    64175// sort two times vectors and an index by first time vector
    65 void sort_mosaic_times (unsigned int *S, unsigned int *E, off_t *I, off_t N) {
    66 
    67 # define SWAPFUNC(A,B){ unsigned int tmp_t; off_t tmp_i; \
     176void sort_image_times (e_time *S, e_time *E, short *C, unsigned int *I, unsigned int *Q, off_t N) {
     177
     178# define SWAPFUNC(A,B){ e_time tmp_t; unsigned int tmp_u; short tmp_c;  \
    68179  tmp_t = S[A]; S[A] = S[B]; S[B] = tmp_t; \
    69180  tmp_t = E[A]; E[A] = E[B]; E[B] = tmp_t; \
    70   tmp_i = I[A]; I[A] = I[B]; I[B] = tmp_i; \
     181  tmp_c = C[A]; C[A] = C[B]; C[B] = tmp_c; \
     182  tmp_u = I[A]; I[A] = I[B]; I[B] = tmp_u; \
     183  tmp_u = Q[A]; Q[A] = Q[B]; Q[B] = tmp_u; \
    71184}
    72185# define COMPARE(A,B)(S[A] < S[B])
    … …  
    79192}
    80193
    81 // use the time to make the match, but use bracketing to make it faster:
    82 
    83 // find mosaic frames (unique time periods & photcode name matches mosaic)
    84 void initMosaics (ImageSubset *image, off_t Nimage) {
    85 
    86   off_t i, Nmos, NMOSAIC;
    87   unsigned int start, stop;
    88 
    89   unsigned int *startMos, *stopMos;
    90   off_t *indexMos;
    91 
    92   Nmosaic = 0;
    93   NMOSAIC = 10;
    94   ALLOCATE (mosaic, Mosaic, NMOSAIC);
    95 
    96   ALLOCATE (startMos, unsigned int, NMOSAIC);
    97   ALLOCATE (stopMos, unsigned int, NMOSAIC);
    98   ALLOCATE (indexMos, off_t, NMOSAIC);
    99 
    100   /* find the mosaic images (coords.ctype = DIS); generate list of unique mosaics */
    101   for (i = 0; i < Nimage; i++) {
    102     if (strcmp(&image[i].coords.ctype[4], "-DIS")) continue;
    103 
    104     /* set image time range */
    105     start = image[i].tzero - MAX(0.01*image[i].trate*image[i].NY, 1);
    106     stop  = image[i].tzero + MAX(1.01*image[i].trate*image[i].NY, 1);
    107 
    108     /* a new mosaic, define ranges */
    109     mosaic[Nmosaic].start = start;
    110     mosaic[Nmosaic].stop  = stop;
    111     mosaic[Nmosaic].coords = image[i].coords;
    112 
    113     startMos[Nmosaic] = start;
    114     stopMos[Nmosaic] = stop;
    115     indexMos[Nmosaic] = Nmosaic;
    116 
    117     Nmosaic ++;
    118     if (Nmosaic == NMOSAIC) {
    119       NMOSAIC += 10;
    120       REALLOCATE (mosaic, Mosaic, NMOSAIC);
    121       REALLOCATE (startMos, unsigned int, NMOSAIC);
    122       REALLOCATE (stopMos, unsigned int, NMOSAIC);
    123       REALLOCATE (indexMos, off_t, NMOSAIC);
    124     }
    125   }
    126 
    127   // sort the index, start, and stop by the start times:
    128   sort_mosaic_times (startMos, stopMos, indexMos, Nmosaic);
    129  
    130   // array to store image->mosaic index
    131   Nmosaic_for_images = Nimage;
    132   ALLOCATE (mosaic_for_images, off_t, Nmosaic_for_images);
    133 
    134   /* now assign the WRP images to these mosaics */
    135   for (i = 0; i < Nimage; i++) {
    136     mosaic_for_images[i] = -1; // default value for no mosaic found
    137 
    138     if (strcmp(&image[i].coords.ctype[4], "-WRP")) continue;
    139 
    140     /* set image time range */
    141     start = image[i].tzero - MAX(0.01*image[i].trate*image[i].NY, 1);
    142     stop  = image[i].tzero + MAX(1.01*image[i].trate*image[i].NY, 1);
    143 
    144     Nmos = getMosaicByTimes (start, stop, startMos, stopMos, indexMos);
    145     if (Nmos == -1) {
    146       fprintf (stderr, "cannot match mosaic for %d (ID = %d)\n", (int) i, (int) image[i].imageID);
    147       continue;
    148     }
    149 
    150     // mosaic corresponding to this image
    151     mosaic_for_images[i] = Nmos;
    152   }
    153 
    154   free (startMos);
    155   free (stopMos);
    156   free (indexMos);
    157   return;
    158 }
    159 
    160 Mosaic *getMosaicForImage (off_t im) {
    161 
    162   off_t mos;
    163 
    164   if (im < 0) abort();
    165   if (im >= Nmosaic_for_images) abort();
    166 
    167   // search for the mosaic that
    168   mos = mosaic_for_images[im];
    169   if (mos < 0) return NULL;
    170 
    171   return &mosaic[mos];
    172 }
    173 
     194# define GET_COLUMN(OUT,NAME,TYPE) \
     195  TYPE *OUT = gfits_get_bintable_column_data (&theader, &ftable, NAME, type, &Nrow, &Ncol); \
     196  myAssert (!strcmp(type, #TYPE), "wrong column type");
     197
     198int ImageValidLoad(char *filename) {
     199
     200  int Ncol;
     201  off_t i;
     202  off_t Nrow;
     203  Header header;
     204  Header theader;
     205  Matrix matrix;
     206  FTable ftable;
     207
     208  FILE *f = fopen (filename, "r");
     209  if (!f) {
     210    fprintf (stderr, "ERROR: cannot open image subset file %s\n", filename);
     211    return FALSE;
     212  }
     213
     214  /* load in PHU segment (ignore) */
     215  if (!gfits_fread_header (f, &header)) {
     216    if (VERBOSE) fprintf (stderr, "can't read image subset header\n");
     217    fclose (f);
     218    return FALSE;
     219  }
     220  if (!gfits_fread_matrix (f, &matrix, &header)) {
     221    if (VERBOSE) fprintf (stderr, "can't read image subset matrix\n");
     222    gfits_free_header (&header);
     223    fclose (f);
     224    return FALSE;
     225  }
     226
     227  ftable.header = &theader;
     228
     229  // load data for this header
     230  if (!gfits_load_header (f, &theader)) {
     231    if (VERBOSE) fprintf (stderr, "can't read table header\n");
     232    fclose (f);
     233    return FALSE;
     234  }
     235
     236  // read the fits table bytes
     237  if (!gfits_fread_ftable_data (f, &ftable, FALSE)) {
     238    if (VERBOSE) fprintf (stderr, "can't read table data\n");
     239    fclose (f);
     240    return FALSE;
     241  }
     242  fclose (f);
     243
     244  char type[16];
     245
     246  GET_COLUMN (Nvalid,     "NVALID",   int);
     247  GET_COLUMN (Ninvalid,   "NINVALID", int);
     248
     249  myAssert (Nrow == Nimage, "wrong size imstats result file\n");
     250
     251  // merge the new arrays into the single array?
     252  for (i = 0; i < Nrow; i++) {
     253    Nimage_valid[i] += Nvalid[i];
     254    Nimage_invalid[i] += Ninvalid[i];
     255  }
     256  fprintf (stderr, "loaded data for %lld images\n", (long long) Nrow);
     257
     258  free (Nvalid);
     259  free (Ninvalid);
     260
     261  return TRUE;
     262}
     263
     264int ImageValidSave(char *filename) {
     265
     266  Header header;
     267  Header theader;
     268  Matrix matrix;
     269  FTable ftable;
     270
     271  gfits_init_header (&header);
     272  header.extend = TRUE;
     273  gfits_create_header (&header);
     274  gfits_create_matrix (&header, &matrix);
     275
     276  gfits_create_table_header (&theader, "BINTABLE", "IMAGE_VALID");
     277
     278  // create the table layout
     279  gfits_define_bintable_column (&theader, "J",   "NVALID",     "tmp",      "tmp", 1.0, FT_BZERO_INT32);
     280  gfits_define_bintable_column (&theader, "J",   "NINVALID",   "tmp",      "tmp", 1.0, FT_BZERO_INT32);
     281
     282  // generate the output array that carries the data
     283  gfits_create_table (&theader, &ftable);
     284
     285  // add the columns to the output array
     286  gfits_set_bintable_column (&theader, &ftable, "NVALID",     Nimage_valid,   Nimage);
     287  gfits_set_bintable_column (&theader, &ftable, "NINVALID",   Nimage_invalid, Nimage);
     288
     289  FILE *f = fopen (filename, "w");
     290  if (!f) {
     291    fprintf (stderr, "ERROR: cannot open image subset file for output %s\n", filename);
     292    return FALSE;
     293  }
     294
     295  gfits_fwrite_header  (f, &header);
     296  gfits_fwrite_matrix  (f, &matrix);
     297  gfits_fwrite_Theader (f, &theader);
     298  gfits_fwrite_table   (f, &ftable);
     299
     300  gfits_free_header (&header);
     301  gfits_free_matrix (&matrix);
     302  gfits_free_header (&theader);
     303  gfits_free_table (&ftable);
     304
     305  fclose (f);
     306  fflush (f);
     307
     308  return TRUE;
     309}

  • branches/eam_branches/ipp-20130509/Ohana/src/uniphot/src/update_catalog_fiximids.c

    r35592 r35593  
    11# include "fiximids.h"
    22
    3 void update_catalog_fiximids (Catalog *catalog, ImageSubset *image, off_t *index, off_t Nimage) {
     3void update_catalog_fiximids (Catalog *catalog) {
    44
    5   float posAngle, pltScale;
    6   off_t i, j, found;
     5  off_t i;
    76
    8   found = 0;   
     7  off_t found = 0;   
     8  off_t Nvalid = 0;   
     9  off_t Ninvalid = 0;   
     10
    911  for (i = 0; i < catalog[0].Nmeasure; i++) {
    1012   
    … …  
    1719    e_time time = catalog[0].measure[i].t;
    1820
    19     off_t imageID = FindImageID (time, photcode);
    20     if (imageID <= 0) continue; // detections with imageID == 0 do not have a valid image (eg, ref photcode)
     21    // skip detections with no valid imageID (eg, ref photcode)
     22
     23    off_t imageID, Seq;
     24    if (!FindImageID (&imageID, &Seq, time, photcode)) continue;
     25    if (imageID <= 0) continue;
    2126
    2227    if (imageID == catalog[0].measure[i].imageID) {
    2328      // existing image ID is valid
    2429      Nvalid ++;
    25       onImage[imageID] ++;
    26       continue;
     30      BumpValidImage (Seq);
     31    } else {
     32      // existing image ID is NOT valid
     33      Ninvalid ++;
     34      catalog[0].measure[i].imageID = imageID;
     35      BumpInvalidImage (Seq);
    2736    }
    28 
    29     // existing image ID is NOT valid
    30     catalog[0].measure[i].imageID = imageID;
    31     Ninvalid ++;
    32     onImage[imageID] ++;
    33     continue;
    3437
    3538    found ++;
    … …  
    3740
    3841  if (found) {
    39     fprintf (stderr, "found "OFF_T_FMT" matches\n", found);
     42    fprintf (stderr, "found "OFF_T_FMT" matches: "OFF_T_FMT" valid, "OFF_T_FMT" invalid\n", found, Nvalid, Ninvalid);
    4043  }
    4144}

  • branches/eam_branches/ipp-20130509/Ohana/src/uniphot/src/update_dvo_fiximids.c

    r35592 r35593  
    66  SkyList *skylist = NULL;
    77  Catalog catalog;
    8   off_t i, maxID, *index;
     8  off_t i;
    99
    1010  // load the current sky table (layout of all SkyRegions)
    … …  
    1212  SkyTableSetFilenames (sky, CATDIR, "cpt");
    1313 
     14  // generate index tables
     15  initImageIndex (image, Nimage);
     16
    1417  if (PARALLEL && !HOST_ID) {
    1518    update_dvo_fiximids_parallel (sky, image, Nimage);
    1619    return TRUE;
    1720  }
    18 
    19   // create an index for the image IDs
    20   maxID = 0;
    21   for (i = 0; i < Nimage; i++) {
    22     maxID = MAX(maxID, image[i].imageID);
    23   }
    24   ALLOCATE (index, off_t, maxID + 1);
    25   for (i = 0; i < maxID + 1; i++) {
    26     index[i] = -1;
    27   }
    28   for (i = 0; i < Nimage; i++) {
    29     if (!image[i].imageID) continue; // images with ID == 0 are virtual
    30     index[image[i].imageID] = i;
    31   }
    32 
    33   // set up the links between images and mosaics
    34   initMosaics (image, Nimage);
    3521
    3622  // determine the populated SkyRegions overlapping the requested area (default depth)
    … …  
    6955    }
    7056
    71     update_catalog_fiximids (&catalog, image, index, Nimage);
     57    update_catalog_fiximids (&catalog);
    7258
    7359    if (!UPDATE) {
    … …  
    8369    dvo_catalog_free (&catalog);
    8470  }
     71
     72  // for the remote clients, we need to dump out the results table (Nvalid, Ninvalid)
     73  if (HOST_ID && IMSTATS_FILE) {
     74    ImageValidSave (IMSTATS_FILE);
     75  }
     76
    8577  return (TRUE);
    8678}     
    … …  
    9082int update_dvo_fiximids_parallel (SkyTable *sky, ImageSubset *image, off_t Nimage) {
    9183
     84  char uniquer[12];
     85  int TIME = time(NULL);
     86  int PID = getpid();
     87  snprintf (uniquer, 12, "%05d.%05d", PID, TIME % 100000);
     88
    9289  // write out the subset table of image information
    9390  char imageFile[512];
    94   snprintf (imageFile, 512, "%s/Images.subset.dat", CATDIR);
     91  snprintf (imageFile, 512, "%s/fiximids.images.%s.dat", CATDIR, uniquer);
    9592
    9693  if (!ImageSubsetSave (imageFile, image, Nimage)) {
    … …  
    112109    table->hosts[i].pathname = tmppath;
    113110
     111    ALLOCATE (table->hosts[i].results, char, 1024);
     112    snprintf (table->hosts[i].results, 1024, "%s/fiximids.results.%s.dat", table->hosts[i].pathname, uniquer);
     113
    114114    char command[1024];
    115115    snprintf (command, 1024, "fiximids_client -hostID %d -catdir %s -hostdir %s -images %s -region %f %f %f %f",
    116               table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, imageFile, 
     116              table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, imageFile,
    117117              UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax);
    118118
    119     fprintf (stderr, "command: %s\n", command);
     119    int saveImstats = VERBOSE_IMSTATS || SUMMARY_IMSTATS;
    120120
    121121    char tmpline[1024];
    122     if (VERBOSE)       { snprintf (tmpline, 1024, "%s -v",                  command); strcpy (command, tmpline); }
    123     if (SINGLE_CPT)    { snprintf (tmpline, 1024, "%s -cpt %s", SINGLE_CPT, command); strcpy (command, tmpline); }
    124     if (UPDATE)        { snprintf (tmpline, 1024, "%s -update",             command); strcpy (command, tmpline); }
     122    if (VERBOSE)     { snprintf (tmpline, 1024, "%s -v",          command);                          strcpy (command, tmpline); }
     123    if (saveImstats) { snprintf (tmpline, 1024, "%s -imstats %s", command, table->hosts[i].results); strcpy (command, tmpline); }
     124    if (SINGLE_CPT)  { snprintf (tmpline, 1024, "%s -cpt %s",     command, SINGLE_CPT);              strcpy (command, tmpline); }
     125    if (UPDATE)      { snprintf (tmpline, 1024, "%s -update",     command);                          strcpy (command, tmpline); }
     126
     127    fprintf (stderr, "command: %s\n", command);
    125128
    126129    if (PARALLEL_MANUAL) continue;
    … …  
    154157  }
    155158
     159  if (VERBOSE_IMSTATS || SUMMARY_IMSTATS) {
     160    for (i = 0; i < table->Nhosts; i++) {
     161      while (!ImageValidLoad (table->hosts[i].results)) {
     162        // failed to get the data from this host.  This can happen for various reasons.  Give the user a chance to try again...
     163        fprintf (stderr, "failed to read data from %s, stopping operations until this can be fixed\n", table->hosts[i].hostname);
     164        fprintf (stderr, "you may run the command manually and send this process the CONT signal\n");
     165        int pid = getpid();
     166        kill (pid, SIGSTOP);
     167        fprintf (stderr, "retrying %s\n", table->hosts[i].results);
     168      }
     169    }
     170  }
     171
    156172  return (TRUE);
    157173}     
Note: See TracChangeset for help on using the changeset viewer.