IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Jul 31, 2012, 4:02:00 PM (14 years ago)
Author:
eugene
Message:

add PS1_V4 schema; add set mean fluxes based on stacks; define the primary image ID (based on distance from center or boundary tree) and set mean flux from stacks based on that info; fix (uniquify) psps ID for stack detections; mean fluxes are in Jy; detection fluxes are instrumental counts per second

Location:
trunk/Ohana/src/addstar
Files:
15 edited
3 copied

Legend:

Unmodified
Added
Removed
  • trunk/Ohana/src/addstar/Makefile

    r34088 r34260  
    1717FULL_LDFLAGS  = -lkapa -ldvo -lFITS -lohana $(BASE_LDFLAGS)
    1818
    19 addstar     : $(BIN)/addstar.$(ARCH)
    20 addstard    : $(BIN)/addstard.$(ARCH)
    21 addstart    : $(BIN)/addstart.$(ARCH)
    22 addstarc    : $(BIN)/addstarc.$(ARCH)
    23 mkacc-2mass : $(BIN)/mkacc-2mass.$(ARCH)
    24 sedstar     : $(BIN)/sedstar.$(ARCH)
    25 load2mass   : $(BIN)/load2mass.$(ARCH)
    26 loadwise    : $(BIN)/loadwise.$(ARCH)
    27 loadsupercos: $(BIN)/loadsupercos.$(ARCH)
    28 gztest      : $(BIN)/gztest.$(ARCH)
    29 mkcmf       : $(BIN)/mkcmf.$(ARCH)
     19addstar      : $(BIN)/addstar.$(ARCH)
     20addstard     : $(BIN)/addstard.$(ARCH)
     21addstart     : $(BIN)/addstart.$(ARCH)
     22addstarc     : $(BIN)/addstarc.$(ARCH)
     23mkacc-2mass  : $(BIN)/mkacc-2mass.$(ARCH)
     24sedstar      : $(BIN)/sedstar.$(ARCH)
     25load2mass    : $(BIN)/load2mass.$(ARCH)
     26loadwise     : $(BIN)/loadwise.$(ARCH)
     27dumpskycells : $(BIN)/dumpskycells.$(ARCH)
     28findskycell  : $(BIN)/findskycell.$(ARCH)
     29loadsupercos : $(BIN)/loadsupercos.$(ARCH)
     30gztest       : $(BIN)/gztest.$(ARCH)
     31mkcmf        : $(BIN)/mkcmf.$(ARCH)
    3032
    3133all: addstar addstar_client sedstar load2mass skycells mkcmf loadwise loadsupercos dumpskycells
    3234
    33 INSTALL = addstar addstar_client sedstar load2mass skycells mkcmf loadwise loadsupercos dumpskycells
     35INSTALL = addstar addstar_client sedstar load2mass skycells mkcmf loadwise loadsupercos dumpskycells findskycell
    3436
    3537# I need to fix the client/server version of addstar now that I have dropped Stars
     
    299301$(SRC)/SetSignals.$(ARCH).o \
    300302
     303FINDSKYCELL = \
     304$(SRC)/findskycell.$(ARCH).o \
     305$(SRC)/Shutdown.$(ARCH).o
     306
    301307$(ADDSTARC)       : $(INC)/addstar.h
    302308$(ADDSTARD)       : $(INC)/addstar.h
     
    306312$(SKYCELLS)       : $(INC)/addstar.h
    307313$(DUMPSKYCELLS)   : $(INC)/addstar.h
     314$(FINDSKYCELL)    : $(INC)/addstar.h
    308315$(LOAD-2MASS)     : $(INC)/addstar.h $(INC)/2mass.h
    309316$(LOAD-WISE)      : $(INC)/addstar.h $(INC)/WISE.h
     
    322329$(BIN)/skycells.$(ARCH)       : $(SKYCELLS)
    323330$(BIN)/dumpskycells.$(ARCH)   : $(DUMPSKYCELLS)
     331$(BIN)/findskycell.$(ARCH)    : $(FINDSKYCELL)
    324332$(BIN)/mkcmf.$(ARCH)          : $(MKCMF)
    325333
  • trunk/Ohana/src/addstar/include/addstar.h

    r33963 r34260  
    297297uint64_t CreatePSPSDetectionID(double tobs, int ccdid, int detID);
    298298uint64_t CreatePSPSObjectID(double ra, double dec);
     299uint64_t CreatePSPSStackDetectionID(int sourceID, int imageID, int detID);
    299300
    300301int altaz (double *alt, double *az, double ha, double dec, double latitude);
  • trunk/Ohana/src/addstar/include/skycells.h

    r33719 r34260  
    1212# include <glob.h>
    1313
    14 enum {SQUARES, TRIANGLES, LOCAL, RINGS};
     14enum {SQUARES, TRIANGLES, LOCAL, RINGS, TAMAS};
    1515enum {TETRAHEDRON, CUBE, OCTOHEDRON, DODECAHEDRON, ICOSAHEDRON};
    1616
     
    9595int          sky_tessellation_squares       PROTO((FITS_DB *db, int level, int Nmax));
    9696int          sky_tessellation_rings         PROTO((FITS_DB *db, int level, int Nmax));
     97int          sky_tessellation_tamas         PROTO((FITS_DB *db, int level, int Nmax));
    9798
    9899int          sky_triangle_to_image          PROTO((Image *image, SkyTriangle *triangle));
     
    104105
    105106SkyRectangle *sky_rectangle_ring            PROTO((float dec, float dDEC, int *nring, char *format));
     107SkyRectangle *sky_rectangle_tamas           PROTO((double *Dec, double dm, double halfa, double halftheta, int *nring, char *format));
    106108
    107109SkyTriangle *sky_divide_triangles           PROTO((SkyTriangle *in, int *ntriangles));
  • trunk/Ohana/src/addstar/src/FilterStars.c

    r29001 r34260  
    8484      stars[N].measure.Map += MTIME - dMs;
    8585    }
     86    if (!isnan(stars[N].measure.Mkron)) {
     87      stars[N].measure.Mkron += MTIME - dMs;
     88    }
     89    if (!isnan(stars[N].measure.FluxPSF)) {
     90      stars[N].measure.FluxPSF /= image[0].exptime;
     91    }
     92    if (!isnan(stars[N].measure.dFluxPSF)) {
     93      stars[N].measure.dFluxPSF /= image[0].exptime;
     94    }
     95    if (!isnan(stars[N].measure.FluxKron)) {
     96      stars[N].measure.FluxKron /= image[0].exptime;
     97    }
     98    if (!isnan(stars[N].measure.dFluxKron)) {
     99      stars[N].measure.dFluxKron /= image[0].exptime;
     100    }
    86101   
    87102    // the external ID is supplied, but do we trust it?
     
    93108      double mjd;
    94109      mjd = ohana_sec_to_mjd (image[0].tzero);
    95       stars[N].measure.extID = CreatePSPSDetectionID(mjd, image[0].ccdnum, stars[N].measure.detID);
     110      int isStack = ((image[0].photcode >= 11000) && (image[0].photcode <= 11400));
     111
     112      if (isStack) {
     113        stars[N].measure.extID = CreatePSPSStackDetectionID(image[0].sourceID, image[0].externID, stars[N].measure.detID);
     114      } else {
     115        stars[N].measure.extID = CreatePSPSDetectionID(mjd, image[0].ccdnum, stars[N].measure.detID);
     116      }
    96117    } else {
    97118      stars[N].measure.extID = 0;
  • trunk/Ohana/src/addstar/src/ReadStarsFITS.c

    r31160 r34260  
    1010  Header theader;
    1111  FTable table;
    12   Stars *stars;
     12  Stars *stars; // Stars contains Average and Measure
    1313 
    1414  if (in_theader == NULL) {
     
    8989    InitStar (&stars[i]);
    9090
    91     stars[i].measure.Xccd      = smpdata[i].X;
    92     stars[i].measure.Yccd      = smpdata[i].Y;
     91    stars[i].measure.Xccd       = smpdata[i].X;
     92    stars[i].measure.Yccd       = smpdata[i].Y;
     93    stars[i].measure.dXccd      = NAN_S_SHORT; // not provided by SMPDATA:
     94    stars[i].measure.dYccd      = NAN_S_SHORT; // not provided by SMPDATA:
     95   
     96    stars[i].measure.posangle   = NAN_S_SHORT; // not provided by SMPDATA:
     97    stars[i].measure.pltscale   = NAN;         // not provided by SMPDATA:
    9398
    9499    if ((smpdata[i].M >= ZeroPt) || isnan(smpdata[i].M)) {
    95       stars[i].measure.M       = NAN;
    96       stars[i].measure.Map     = NAN;
    97     } else {
    98       stars[i].measure.M       = smpdata[i].M;
    99       stars[i].measure.Map     = smpdata[i].M;
    100     }
    101 
    102     stars[i].measure.dM        = smpdata[i].dM*0.001;
    103 
     100      stars[i].measure.M        = NAN;
     101      stars[i].measure.Map      = NAN;
     102      stars[i].measure.FluxPSF  = NAN;
     103      stars[i].measure.dFluxPSF = NAN;
     104    } else {
     105      stars[i].measure.M        = smpdata[i].M;
     106      stars[i].measure.Map      = smpdata[i].M;
     107      stars[i].measure.FluxPSF  = pow(10.0, -0.4*smpdata[i].M);
     108      stars[i].measure.dFluxPSF = stars[i].measure.FluxPSF * smpdata[i].dM;
     109    }
     110    stars[i].measure.dM         = smpdata[i].dM*0.001;
     111    stars[i].measure.dMcal      = NAN; // not provided by SMPDATA:
     112
     113    stars[i].measure.Mkron      = NAN; // not provided by SMPDATA:
     114    stars[i].measure.dMkron     = NAN; // not provided by SMPDATA:
     115    stars[i].measure.FluxKron   = NAN; // not provided by SMPDATA:
     116    stars[i].measure.dFluxKron  = NAN; // not provided by SMPDATA:
     117
     118    stars[i].measure.Sky        = NAN; // not provided by SMPDATA:
     119    stars[i].measure.dSky       = NAN; // not provided by SMPDATA:
     120
     121    stars[i].measure.psfChisq   = NAN;       // not provided by SMPDATA:
     122    stars[i].measure.psfQual    = NAN;       // not provided by SMPDATA:
     123    stars[i].measure.psfNdof    = NAN_S_INT; // not provided by SMPDATA:
     124    stars[i].measure.psfNpix    = NAN_S_INT; // not provided by SMPDATA:
     125    stars[i].measure.crNsigma   = NAN;       // not provided by SMPDATA:
     126    stars[i].measure.extNsigma  = NAN;       // not provided by SMPDATA:
     127
     128    stars[i].measure.FWx        = ToShortPixels (smpdata[i].fx);
     129    stars[i].measure.FWy        = ToShortPixels (smpdata[i].fy);
     130    stars[i].measure.theta      = ToShortDegrees (smpdata[i].df);
     131
     132    stars[i].measure.Mxx        = NAN_S_SHORT; // not provided by SMPDATA:
     133    stars[i].measure.Mxy        = NAN_S_SHORT; // not provided by SMPDATA:
     134    stars[i].measure.Myy        = NAN_S_SHORT; // not provided by SMPDATA:
     135                       
    104136    // the dophot type information gets pushed into the upper 2 bytes of photFlags
    105     stars[i].measure.photFlags = (smpdata[i].dophot << 16);
    106 
    107     stars[i].measure.FWx       = ToShortPixels (smpdata[i].fx);
    108     stars[i].measure.FWy       = ToShortPixels (smpdata[i].fy);
    109     stars[i].measure.theta     = ToShortDegrees (smpdata[i].df);
     137    stars[i].measure.photFlags  = (smpdata[i].dophot << 16);
    110138  }   
    111139  *nstars = Nstars;
     
    131159  for (i = 0; i < Nstars; i++) {
    132160    InitStar (&stars[i]);
    133     stars[i].measure.Xccd     = ps1data[i].X;
    134     stars[i].measure.Yccd     = ps1data[i].Y;
    135 
    136     stars[i].measure.dXccd    = ToShortPixels(ps1data[i].dX);
    137     stars[i].measure.dYccd    = ToShortPixels(ps1data[i].dY);
     161    stars[i].measure.Xccd       = ps1data[i].X;
     162    stars[i].measure.Yccd       = ps1data[i].Y;
     163    stars[i].measure.dXccd      = ToShortPixels(ps1data[i].dX);
     164    stars[i].measure.dYccd      = ToShortPixels(ps1data[i].dY);
    138165   
     166    stars[i].measure.posangle   = NAN_S_SHORT; // not provided by PS1_DEV_0:
     167    stars[i].measure.pltscale   = NAN;         // not provided by PS1_DEV_0:
     168
    139169    if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    140       stars[i].measure.M      = NAN;
    141     } else {
    142       stars[i].measure.M      = ps1data[i].M + ZeroPt;
    143     }
    144     stars[i].measure.Map      = NAN;
    145     stars[i].measure.dM       = ps1data[i].dM;
    146     stars[i].measure.Sky      = ps1data[i].sky;
    147     stars[i].measure.dSky     = ps1data[i].dSky;
    148 
    149     stars[i].measure.FWx      = ToShortPixels(ps1data[i].fx);
    150     stars[i].measure.FWy      = ToShortPixels(ps1data[i].fy);
    151     stars[i].measure.theta    = ToShortDegrees(ps1data[i].df);
    152 
    153     stars[i].measure.psfChisq = ps1data[i].psfChisq;
    154     stars[i].measure.psfQual  = ps1data[i].psfQual;
    155 
    156     stars[i].measure.detID    = ps1data[i].detID;
     170      stars[i].measure.M        = NAN;
     171      stars[i].measure.FluxPSF  = NAN;
     172      stars[i].measure.dFluxPSF = NAN;
     173    } else {
     174      stars[i].measure.M        = ps1data[i].M + ZeroPt;
     175      stars[i].measure.FluxPSF  = pow(10.0, -0.4*ps1data[i].M);
     176      stars[i].measure.dFluxPSF = stars[i].measure.FluxPSF * ps1data[i].dM;
     177    }
     178    stars[i].measure.dM         = ps1data[i].dM;
     179    stars[i].measure.dMcal      = NAN; // not provided by PS1_DEV_0:
     180    stars[i].measure.Map        = NAN; // not provided by PS1_DEV_0:
     181
     182    stars[i].measure.Mkron      = NAN; // not provided by PS1_DEV_0:
     183    stars[i].measure.dMkron     = NAN; // not provided by PS1_DEV_0:
     184    stars[i].measure.FluxKron   = NAN; // not provided by PS1_DEV_0:
     185    stars[i].measure.dFluxKron  = NAN; // not provided by PS1_DEV_0:
     186
     187    stars[i].measure.Sky        = ps1data[i].sky;
     188    stars[i].measure.dSky       = ps1data[i].dSky;
     189
     190    stars[i].measure.psfChisq   = ps1data[i].psfChisq;
     191    stars[i].measure.psfQual    = ps1data[i].psfQual;
     192    stars[i].measure.psfNdof    = NAN_S_INT; // not provided by PS1_DEV_0:
     193    stars[i].measure.psfNpix    = NAN_S_INT; // not provided by PS1_DEV_0:
     194    stars[i].measure.crNsigma   = NAN;       // not provided by PS1_DEV_0:
     195    stars[i].measure.extNsigma  = NAN;        // not provided by PS1_DEV_0:
     196
     197    stars[i].measure.FWx        = ToShortPixels(ps1data[i].fx);
     198    stars[i].measure.FWy        = ToShortPixels(ps1data[i].fy);
     199    stars[i].measure.theta      = ToShortDegrees(ps1data[i].df);
     200
     201    stars[i].measure.Mxx        = NAN_S_SHORT; // not provided by PS1_DEV_0:
     202    stars[i].measure.Mxy        = NAN_S_SHORT; // not provided by PS1_DEV_0:
     203    stars[i].measure.Myy        = NAN_S_SHORT; // not provided by PS1_DEV_0:
     204                       
     205    stars[i].measure.photFlags  = 0; // not provided by PS1_DEV_0:
     206
     207    stars[i].measure.detID      = ps1data[i].detID;
    157208  }   
    158209  *nstars = Nstars;
     
    182233    stars[i].measure.Xccd       = ps1data[i].X;
    183234    stars[i].measure.Yccd       = ps1data[i].Y;
    184 
    185235    stars[i].measure.dXccd      = ToShortPixels(ps1data[i].dX);
    186236    stars[i].measure.dYccd      = ToShortPixels(ps1data[i].dY);
    187237
     238    stars[i].measure.posangle   = NAN_S_SHORT; // not provided by PS1_DEV_1:
     239    stars[i].measure.pltscale   = NAN;         // not provided by PS1_DEV_1:
     240
    188241    if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    189         stars[i].measure.M      = NAN;
    190     } else {
    191         stars[i].measure.M      = ps1data[i].M + ZeroPt;
    192     }
    193     stars[i].measure.Map        = NAN;
     242      stars[i].measure.M      = NAN;
     243      stars[i].measure.FluxPSF  = NAN;
     244      stars[i].measure.dFluxPSF = NAN;
     245    } else {
     246      stars[i].measure.M      = ps1data[i].M + ZeroPt;
     247      stars[i].measure.FluxPSF  = pow(10.0, -0.4*ps1data[i].M);
     248      stars[i].measure.dFluxPSF = stars[i].measure.FluxPSF * ps1data[i].dM;
     249    }
    194250    stars[i].measure.dM         = ps1data[i].dM;
     251    stars[i].measure.dMcal      = NAN; // not provided by PS1_DEV_1:
     252    stars[i].measure.Map        = NAN; // not provided by PS1_DEV_1:
     253
     254    stars[i].measure.Mkron      = NAN; // not provided by PS1_DEV_1:
     255    stars[i].measure.dMkron     = NAN; // not provided by PS1_DEV_1:
     256    stars[i].measure.FluxKron   = NAN; // not provided by PS1_DEV_1:
     257    stars[i].measure.dFluxKron  = NAN; // not provided by PS1_DEV_1:
     258
    195259    stars[i].measure.Sky        = ps1data[i].sky;
    196260    stars[i].measure.dSky       = ps1data[i].dSky;
     261
     262    stars[i].measure.psfChisq   = ps1data[i].psfChisq;
     263    stars[i].measure.psfQual    = ps1data[i].psfQual;
     264    stars[i].measure.psfNdof    = NAN_S_INT; // not provided by PS1_DEV_1:
     265    stars[i].measure.psfNpix    = NAN_S_INT; // not provided by PS1_DEV_1:
     266    stars[i].measure.crNsigma   = ps1data[i].crNsigma;
     267    stars[i].measure.extNsigma  = ps1data[i].extNsigma;
    197268
    198269    stars[i].measure.FWx        = ToShortPixels(ps1data[i].fx);
     
    200271    stars[i].measure.theta      = ToShortDegrees(ps1data[i].df);
    201272
    202     stars[i].measure.psfChisq   = ps1data[i].psfChisq;
    203     stars[i].measure.psfQual    = ps1data[i].psfQual;
    204     stars[i].measure.crNsigma   = ps1data[i].crNsigma;
    205     stars[i].measure.extNsigma  = ps1data[i].extNsigma;
    206 
    207     stars[i].measure.detID      = ps1data[i].detID;
     273    stars[i].measure.Mxx        = NAN_S_SHORT; // not provided by PS1_DEV_1:
     274    stars[i].measure.Mxy        = NAN_S_SHORT; // not provided by PS1_DEV_1:
     275    stars[i].measure.Myy        = NAN_S_SHORT; // not provided by PS1_DEV_1:
     276                       
    208277    stars[i].measure.photFlags  = ps1data[i].flags;
     278
     279    // this is may optionally be replaced by the internal sequence (see FilterStars.c)
     280    stars[i].measure.detID      = ps1data[i].detID;
    209281  }   
    210282  *nstars = Nstars;
     
    229301
    230302  if (table[0].header[0].Naxis[0] == 136) {
    231       stars = Convert_PS1_V1_Alt (table, nstars);
    232       return (stars);
     303    stars = Convert_PS1_V1_Alt (table, nstars);
     304    return (stars);
    233305  }
    234306
     
    252324
    253325    if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    254         stars[i].measure.M      = NAN;
    255     } else {
    256         stars[i].measure.M      = ps1data[i].M + ZeroPt;
     326      stars[i].measure.M      = NAN;
     327      stars[i].measure.FluxPSF  = NAN;
     328      stars[i].measure.dFluxPSF = NAN;
     329    } else {
     330      stars[i].measure.M      = ps1data[i].M + ZeroPt;
     331      stars[i].measure.FluxPSF  = pow(10.0, -0.4*ps1data[i].M);
     332      stars[i].measure.dFluxPSF = stars[i].measure.FluxPSF * ps1data[i].dM;
    257333    }
    258334    stars[i].measure.dM         = ps1data[i].dM;
    259335    stars[i].measure.dMcal      = ps1data[i].dMcal;
    260336    stars[i].measure.Map        = ps1data[i].Map + ZeroPt;
    261                        
     337                       
     338    stars[i].measure.Mkron      = NAN; // not provided by PS1_V1:
     339    stars[i].measure.dMkron     = NAN; // not provided by PS1_V1:
     340    stars[i].measure.FluxKron   = NAN; // not provided by PS1_V1:
     341    stars[i].measure.dFluxKron  = NAN; // not provided by PS1_V1:
     342
    262343    stars[i].measure.Sky        = ps1data[i].sky;
    263344    stars[i].measure.dSky       = ps1data[i].dSky;
    264                        
     345                       
    265346    stars[i].measure.psfChisq   = ps1data[i].psfChisq;
    266347    stars[i].measure.psfQual    = ps1data[i].psfQual;
     
    277358    stars[i].measure.Mxy        = ToShortPixels(ps1data[i].Mxy);
    278359    stars[i].measure.Myy        = ToShortPixels(ps1data[i].Myy);
    279                        
     360                       
    280361    stars[i].measure.photFlags  = ps1data[i].flags;
    281362
     
    328409
    329410    if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    330         stars[i].measure.M      = NAN;
    331     } else {
    332         stars[i].measure.M      = ps1data[i].M + ZeroPt;
     411      stars[i].measure.M      = NAN;
     412      stars[i].measure.FluxPSF  = NAN;
     413      stars[i].measure.dFluxPSF = NAN;
     414    } else {
     415      stars[i].measure.M        = ps1data[i].M + ZeroPt;
     416      stars[i].measure.FluxPSF  = pow(10.0, -0.4*ps1data[i].M);
     417      stars[i].measure.dFluxPSF = stars[i].measure.FluxPSF * ps1data[i].dM;
    333418    }
    334419    stars[i].measure.dM         = ps1data[i].dM;
    335420    stars[i].measure.dMcal      = ps1data[i].dMcal;
    336421    stars[i].measure.Map        = ps1data[i].Map + ZeroPt;
    337                        
     422                       
     423    stars[i].measure.Mkron      = NAN; // not provided by PS1_V1_Alt:
     424    stars[i].measure.dMkron     = NAN; // not provided by PS1_V1_Alt:
     425    stars[i].measure.FluxKron   = NAN; // not provided by PS1_V1_Alt:
     426    stars[i].measure.dFluxKron  = NAN; // not provided by PS1_V1_Alt:
     427
    338428    stars[i].measure.Sky        = ps1data[i].sky;
    339429    stars[i].measure.dSky       = ps1data[i].dSky;
    340                        
     430                       
    341431    stars[i].measure.psfChisq   = ps1data[i].psfChisq;
    342432    stars[i].measure.psfQual    = ps1data[i].psfQual;
     
    353443    stars[i].measure.Mxy        = ToShortPixels(ps1data[i].Mxy);
    354444    stars[i].measure.Myy        = ToShortPixels(ps1data[i].Myy);
    355                        
     445                       
    356446    stars[i].measure.photFlags  = ps1data[i].flags;
    357447
     
    396486
    397487    if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    398         stars[i].measure.M      = NAN;
    399     } else {
    400         stars[i].measure.M      = ps1data[i].M + ZeroPt;
     488      stars[i].measure.M        = NAN;
     489      stars[i].measure.FluxPSF  = NAN;
     490      stars[i].measure.dFluxPSF = NAN;
     491    } else {
     492      stars[i].measure.M      = ps1data[i].M + ZeroPt;
     493      stars[i].measure.FluxPSF  = pow(10.0, -0.4*ps1data[i].M);
     494      stars[i].measure.dFluxPSF = stars[i].measure.FluxPSF * ps1data[i].dM;
    401495    }
    402496    stars[i].measure.dM         = ps1data[i].dM;
    403497    stars[i].measure.dMcal      = ps1data[i].dMcal;
    404498    stars[i].measure.Map        = ps1data[i].Map + ZeroPt;
    405                        
     499                       
     500    stars[i].measure.Mkron      = NAN; // not provided by PS1_V2:
     501    stars[i].measure.dMkron     = NAN; // not provided by PS1_V2:
     502    stars[i].measure.FluxKron   = NAN; // not provided by PS1_V2:
     503    stars[i].measure.dFluxKron  = NAN; // not provided by PS1_V2:
     504
    406505    stars[i].measure.Sky        = ps1data[i].sky;
    407506    stars[i].measure.dSky       = ps1data[i].dSky;
    408                        
     507                       
    409508    stars[i].measure.psfChisq   = ps1data[i].psfChisq;
    410509    stars[i].measure.psfQual    = ps1data[i].psfQual;
     
    421520    stars[i].measure.Mxy        = ToShortPixels(ps1data[i].Mxy);
    422521    stars[i].measure.Myy        = ToShortPixels(ps1data[i].Myy);
    423                        
     522                       
    424523    stars[i].measure.photFlags  = ps1data[i].flags;
    425524
     
    464563
    465564    if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    466         stars[i].measure.M      = NAN;
    467     } else {
    468         stars[i].measure.M      = ps1data[i].M + ZeroPt;
     565      stars[i].measure.M      = NAN;
     566    } else {
     567      stars[i].measure.M      = ps1data[i].M + ZeroPt;
    469568    }
    470569    stars[i].measure.dM         = ps1data[i].dM;
    471570    stars[i].measure.dMcal      = ps1data[i].dMcal;
    472571    stars[i].measure.Map        = ps1data[i].Map + ZeroPt;
    473                        
     572                       
     573    stars[i].measure.Mkron      = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN;
     574    stars[i].measure.dMkron     = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN;
     575                       
     576    // these fluxes are converted from counts to counts/sec in FilterStars.c
     577    stars[i].measure.FluxPSF    = ps1data[i].Flux;
     578    stars[i].measure.dFluxPSF   = ps1data[i].dFlux;
     579    stars[i].measure.FluxKron   = ps1data[i].kronFlux;
     580    stars[i].measure.dFluxKron  = ps1data[i].kronFluxErr;
     581
    474582    stars[i].measure.Sky        = ps1data[i].sky;
    475583    stars[i].measure.dSky       = ps1data[i].dSky;
    476                        
     584                       
    477585    stars[i].measure.psfChisq   = ps1data[i].psfChisq;
    478586    stars[i].measure.psfQual    = ps1data[i].psfQual;
     
    489597    stars[i].measure.Mxy        = ToShortPixels(ps1data[i].Mxy);
    490598    stars[i].measure.Myy        = ToShortPixels(ps1data[i].Myy);
    491                        
     599                       
    492600    stars[i].measure.photFlags  = ps1data[i].flags;
    493601
     
    496604
    497605    // the Average fields and the following Measure fields are set in FilterStars after
    498     // the image metadata is in hand:  dR, dD, Mcal, dt, airmass, az, t, imageID, extID,
    499     // averef is set in find_matches, dbFlags is zero on ingest.
     606    // the image metadata is in hand:  dR, dD, Mcal, dt, airmass, az, t, imageID, extID.
     607
     608    // averef is set in find_matches
     609
     610    // dbFlags is zero on ingest.
    500611
    501612    // the following fields are currently not being set anywhere: t_msec
     
    514625
    515626  if (table[0].header[0].Naxis[0] == 196) {
    516       stars = Convert_PS1_SV1_Alt (table, nstars);
    517       return (stars);
     627    stars = Convert_PS1_SV1_Alt (table, nstars);
     628    return (stars);
    518629  }
    519630
     
    537648
    538649    if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    539         stars[i].measure.M      = NAN;
    540     } else {
    541         stars[i].measure.M      = ps1data[i].M + ZeroPt;
     650      stars[i].measure.M      = NAN;
     651    } else {
     652      stars[i].measure.M      = ps1data[i].M + ZeroPt;
    542653    }
    543654    stars[i].measure.dM         = ps1data[i].dM;
    544655    stars[i].measure.dMcal      = ps1data[i].dMcal;
    545656    stars[i].measure.Map        = ps1data[i].Map + ZeroPt;
    546                        
     657                       
     658    stars[i].measure.Mkron      = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN;
     659    stars[i].measure.dMkron     = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN;
     660
     661    // these fluxes are converted from counts to counts/sec in FilterStars.c
     662    stars[i].measure.FluxPSF    = ps1data[i].Flux;
     663    stars[i].measure.dFluxPSF   = ps1data[i].dFlux;
     664    stars[i].measure.FluxKron   = ps1data[i].kronFlux;
     665    stars[i].measure.dFluxKron  = ps1data[i].kronFluxErr;
     666
    547667    stars[i].measure.Sky        = ps1data[i].sky;
    548668    stars[i].measure.dSky       = ps1data[i].dSky;
    549                        
     669                       
    550670    stars[i].measure.psfChisq   = ps1data[i].psfChisq;
    551671    stars[i].measure.psfQual    = ps1data[i].psfQual;
     
    562682    stars[i].measure.Mxy        = ToShortPixels(ps1data[i].Mxy);
    563683    stars[i].measure.Myy        = ToShortPixels(ps1data[i].Myy);
    564                        
     684                       
    565685    stars[i].measure.photFlags  = ps1data[i].flags;
    566686
     
    607727
    608728    if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    609         stars[i].measure.M      = NAN;
    610     } else {
    611         stars[i].measure.M      = ps1data[i].M + ZeroPt;
     729      stars[i].measure.M      = NAN;
     730    } else {
     731      stars[i].measure.M      = ps1data[i].M + ZeroPt;
    612732    }
    613733    stars[i].measure.dM         = ps1data[i].dM;
    614734    stars[i].measure.dMcal      = ps1data[i].dMcal;
    615735    stars[i].measure.Map        = ps1data[i].Map + ZeroPt;
    616                        
     736                       
     737    stars[i].measure.Mkron      = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN;
     738    stars[i].measure.dMkron     = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN;
     739
     740    // these fluxes are converted from counts to counts/sec in FilterStars.c
     741    stars[i].measure.FluxPSF    = ps1data[i].Flux;
     742    stars[i].measure.dFluxPSF   = ps1data[i].dFlux;
     743    stars[i].measure.FluxKron   = ps1data[i].kronFlux;
     744    stars[i].measure.dFluxKron  = ps1data[i].kronFluxErr;
     745
    617746    stars[i].measure.Sky        = ps1data[i].sky;
    618747    stars[i].measure.dSky       = ps1data[i].dSky;
    619                        
     748                       
    620749    stars[i].measure.psfChisq   = ps1data[i].psfChisq;
    621750    stars[i].measure.psfQual    = ps1data[i].psfQual;
     
    632761    stars[i].measure.Mxy        = ToShortPixels(ps1data[i].Mxy);
    633762    stars[i].measure.Myy        = ToShortPixels(ps1data[i].Myy);
    634                        
     763                       
    635764    stars[i].measure.photFlags  = ps1data[i].flags;
    636765
  • trunk/Ohana/src/addstar/src/StarOps.c

    r30613 r34260  
    33int InitStar (Stars *star) {
    44
    5     memset (&star[0].average, 0, sizeof(Average));
    6     memset (&star[0].measure, 0, sizeof(Measure));
     5
     6    dvo_measure_init (&star[0].measure);
     7    dvo_average_init (&star[0].average);
    78    star[0].found = -1; // found == -1 -> not yet found (use enums?)
    89
  • trunk/Ohana/src/addstar/src/args_skycells.c

    r31239 r34260  
    3737    if (!strcasecmp (argv[N], "rings")) {
    3838      MODE = RINGS;
     39    }
     40    if (!strcasecmp (argv[N], "tamas")) {
     41      MODE = TAMAS;
    3942    }
    4043    remove_argument (N, &argc, argv);
     
    179182    } 
    180183    remove_argument (N, &argc, argv);
     184  }
     185  if (MODE == TAMAS) {
     186    CELLSIZE = 3.955;
     187    if ((N = get_argument (argc, argv, "-cellsize"))) {
     188      remove_argument (N, &argc, argv);
     189      CELLSIZE = strtod (argv[N], &ptr);
     190      if ((*ptr != 0) || (CELLSIZE < 0.0)) {
     191        fprintf (stderr, "-cellsize requires a floating-point argument\n");
     192        help ();
     193      } 
     194      remove_argument (N, &argc, argv);
     195    }
    181196  }
    182197
  • trunk/Ohana/src/addstar/src/find_matches.c

    r33963 r34260  
    250250    if (!IN_REGION (stars[i].average.R, stars[i].average.D)) continue;
    251251
     252    dvo_average_init (&catalog[0].average[Nave]);
    252253    catalog[0].average[Nave].R             = stars[i].average.R;
    253254    catalog[0].average[Nave].D             = stars[i].average.D;
    254     catalog[0].average[Nave].dR            = 0;
    255     catalog[0].average[Nave].dD            = 0;
    256255
    257256    catalog[0].average[Nave].Nmeasure      = NSTAR_GROUP;
    258     catalog[0].average[Nave].Nmissing      = 0;
    259     catalog[0].average[Nave].Nextend       = 0;
    260 
    261257    catalog[0].average[Nave].measureOffset = Nmeas;
    262     catalog[0].average[Nave].missingOffset = -1;
    263     catalog[0].average[Nave].extendOffset  = -1;
    264 
    265     catalog[0].average[Nave].uR            = 0;
    266     catalog[0].average[Nave].uD            = 0;
    267     catalog[0].average[Nave].duR           = 0;
    268     catalog[0].average[Nave].duD           = 0;
    269     catalog[0].average[Nave].P             = 0;
    270     catalog[0].average[Nave].dP            = 0;
    271 
    272     catalog[0].average[Nave].Xp            = 0;
    273     catalog[0].average[Nave].ChiSqAve      = 0.0;
    274     catalog[0].average[Nave].ChiSqPM       = 0.0;
    275     catalog[0].average[Nave].ChiSqPar      = 0.0;
    276     catalog[0].average[Nave].Tmean         = 0;
    277     catalog[0].average[Nave].Trange        = 0;
    278     catalog[0].average[Nave].Npos          = 0;
    279 
    280258    catalog[0].average[Nave].objID         = objID;
    281259    catalog[0].average[Nave].catID         = catID;
    282     catalog[0].average[Nave].flags         = 0;
     260
    283261    if (PSPS_ID) {
    284         catalog[0].average[Nave].extID = CreatePSPSObjectID(catalog[0].average[Nave].R,
    285                                                             catalog[0].average[Nave].D);
    286     } else {
    287         catalog[0].average[Nave].extID         = 0;
     262        catalog[0].average[Nave].extID = CreatePSPSObjectID(catalog[0].average[Nave].R, catalog[0].average[Nave].D);
    288263    }
    289264
     
    291266
    292267    for (j = 0; j < Nsecfilt; j++) {
    293       catalog[0].secfilt[Nave*Nsecfilt+j].M           = NAN;
    294       catalog[0].secfilt[Nave*Nsecfilt+j].Map         = NAN;
    295       catalog[0].secfilt[Nave*Nsecfilt+j].dM          = NAN;
    296       catalog[0].secfilt[Nave*Nsecfilt+j].Mstdev      = NAN_S_SHORT;
    297       catalog[0].secfilt[Nave*Nsecfilt+j].Xm          = NAN_S_SHORT;
    298       catalog[0].secfilt[Nave*Nsecfilt+j].M_20        = NAN_S_SHORT;
    299       catalog[0].secfilt[Nave*Nsecfilt+j].M_80        = NAN_S_SHORT;
    300       catalog[0].secfilt[Nave*Nsecfilt+j].Ncode       = 0;
    301       catalog[0].secfilt[Nave*Nsecfilt+j].Nused       = 0;
    302       catalog[0].secfilt[Nave*Nsecfilt+j].ubercalDist = 1000;
    303       catalog[0].secfilt[Nave*Nsecfilt+j].flags       = 0;
     268      dvo_secfilt_init (&catalog[0].secfilt[Nave*Nsecfilt+j]);
    304269    }
    305270
  • trunk/Ohana/src/addstar/src/find_matches_closest.c

    r33963 r34260  
    252252    if (!IN_REGION (stars[i].average.R, stars[i].average.D)) continue;
    253253
     254    dvo_average_init (&catalog[0].average[Nave]);
    254255    catalog[0].average[Nave].R             = stars[i].average.R;
    255256    catalog[0].average[Nave].D             = stars[i].average.D;
    256     catalog[0].average[Nave].dR            = 0;
    257     catalog[0].average[Nave].dD            = 0;
    258257
    259258    catalog[0].average[Nave].Nmeasure      = NSTAR_GROUP;
    260     catalog[0].average[Nave].Nmissing      = 0;
    261     catalog[0].average[Nave].Nextend       = 0;
    262 
    263259    catalog[0].average[Nave].measureOffset = Nmeas;
    264     catalog[0].average[Nave].missingOffset = -1;
    265     catalog[0].average[Nave].extendOffset  = -1;
    266 
    267     catalog[0].average[Nave].uR            = 0;
    268     catalog[0].average[Nave].uD            = 0;
    269     catalog[0].average[Nave].duR           = 0;
    270     catalog[0].average[Nave].duD           = 0;
    271     catalog[0].average[Nave].P             = 0;
    272     catalog[0].average[Nave].dP            = 0;
    273 
    274     catalog[0].average[Nave].Xp            = 0;
    275     catalog[0].average[Nave].ChiSqAve      = 0.0;
    276     catalog[0].average[Nave].ChiSqPM       = 0.0;
    277     catalog[0].average[Nave].ChiSqPar      = 0.0;
    278     catalog[0].average[Nave].Tmean         = 0;
    279     catalog[0].average[Nave].Trange        = 0;
    280     catalog[0].average[Nave].Npos          = 0;
    281 
    282260    catalog[0].average[Nave].objID         = objID;
    283261    catalog[0].average[Nave].catID         = catID;
    284     catalog[0].average[Nave].flags         = 0;
     262
    285263    if (PSPS_ID) {
    286264        catalog[0].average[Nave].extID = CreatePSPSObjectID(catalog[0].average[Nave].R, catalog[0].average[Nave].D);
    287     } else {
    288         catalog[0].average[Nave].extID = 0;
    289265    }
    290266
     
    292268
    293269    for (j = 0; j < Nsecfilt; j++) {
    294       catalog[0].secfilt[Nave*Nsecfilt+j].M           = NAN;
    295       catalog[0].secfilt[Nave*Nsecfilt+j].Map         = NAN;
    296       catalog[0].secfilt[Nave*Nsecfilt+j].dM          = NAN;
    297       catalog[0].secfilt[Nave*Nsecfilt+j].Mstdev      = NAN_S_SHORT;
    298       catalog[0].secfilt[Nave*Nsecfilt+j].Xm          = NAN_S_SHORT;
    299       catalog[0].secfilt[Nave*Nsecfilt+j].M_20        = NAN_S_SHORT;
    300       catalog[0].secfilt[Nave*Nsecfilt+j].M_80        = NAN_S_SHORT;
    301       catalog[0].secfilt[Nave*Nsecfilt+j].Ncode       = 0;
    302       catalog[0].secfilt[Nave*Nsecfilt+j].Nused       = 0;
    303       catalog[0].secfilt[Nave*Nsecfilt+j].ubercalDist = 1000;
    304       catalog[0].secfilt[Nave*Nsecfilt+j].flags       = 0;
     270      dvo_secfilt_init (&catalog[0].secfilt[Nave*Nsecfilt+j]);
    305271    }
    306272
  • trunk/Ohana/src/addstar/src/find_matches_closest_refstars.c

    r28241 r34260  
    255255    if (!IN_REGION (stars[N][0].average.R, stars[N][0].average.D)) continue;
    256256
     257    dvo_average_init (&catalog[0].average[Nave]);
    257258    catalog[0].average[Nave].R             = stars[N][0].average.R;
    258259    catalog[0].average[Nave].D             = stars[N][0].average.D;
    259260
    260261    catalog[0].average[Nave].Nmeasure      = NREFSTAR_GROUP;
    261     catalog[0].average[Nave].Nmissing      = 0;
    262     catalog[0].average[Nave].Nextend       = 0;
    263 
    264262    catalog[0].average[Nave].measureOffset = Nmeas;
    265     catalog[0].average[Nave].missingOffset = -1;
    266     catalog[0].average[Nave].extendOffset  = -1;
     263    catalog[0].average[Nave].objID         = objID;
     264    catalog[0].average[Nave].catID         = catID;
     265
     266    if (PSPS_ID) {
     267        catalog[0].average[Nave].extID = CreatePSPSObjectID(catalog[0].average[Nave].R, catalog[0].average[Nave].D);
     268    }
    267269
    268270    if (ACCEPT_MOTION) {
     
    275277      catalog[0].average[Nave].P           = stars[N][0].average.P;
    276278      catalog[0].average[Nave].dP          = stars[N][0].average.dP;
    277     } else {
    278       catalog[0].average[Nave].dR          = 0;
    279       catalog[0].average[Nave].dD          = 0;
    280       catalog[0].average[Nave].uR          = 0;
    281       catalog[0].average[Nave].uD          = 0;
    282       catalog[0].average[Nave].duR         = 0;
    283       catalog[0].average[Nave].duD         = 0;
    284       catalog[0].average[Nave].P           = 0;
    285       catalog[0].average[Nave].dP          = 0;
    286       catalog[0].average[Nave].Xp          = 0;
    287     }
    288 
    289     catalog[0].average[Nave].Xp            = 0;
    290     catalog[0].average[Nave].ChiSqAve      = 0.0;
    291     catalog[0].average[Nave].ChiSqPM       = 0.0;
    292     catalog[0].average[Nave].ChiSqPar      = 0.0;
    293     catalog[0].average[Nave].Tmean         = 0;
    294     catalog[0].average[Nave].Trange        = 0;
    295     catalog[0].average[Nave].Npos          = 0;
    296 
    297     catalog[0].average[Nave].objID         = objID;
    298     catalog[0].average[Nave].catID         = catID;
    299     catalog[0].average[Nave].flags         = 0;
    300     if (PSPS_ID) {
    301         catalog[0].average[Nave].extID = CreatePSPSObjectID(catalog[0].average[Nave].R,
    302                                                             catalog[0].average[Nave].D);
    303     } else {
    304         catalog[0].average[Nave].extID         = 0;
    305279    }
    306280
     
    308282
    309283    for (j = 0; j < Nsecfilt; j++) {
    310       catalog[0].secfilt[Nave*Nsecfilt+j].M     = NAN;
    311       catalog[0].secfilt[Nave*Nsecfilt+j].dM    = NAN;
    312       catalog[0].secfilt[Nave*Nsecfilt+j].Xm    = NAN_S_SHORT;
    313       catalog[0].secfilt[Nave*Nsecfilt+j].M_20  = NAN_S_SHORT;
    314       catalog[0].secfilt[Nave*Nsecfilt+j].M_80  = NAN_S_SHORT;
    315       catalog[0].secfilt[Nave*Nsecfilt+j].Ncode = 0;
    316       catalog[0].secfilt[Nave*Nsecfilt+j].Nused = 0;
     284      dvo_secfilt_init (&catalog[0].secfilt[Nave*Nsecfilt+j]);
    317285    }
    318286
  • trunk/Ohana/src/addstar/src/find_matches_refstars.c

    r33653 r34260  
    227227    if (!IN_REGION (stars[N][0].average.R, stars[N][0].average.D)) continue;
    228228
     229    dvo_average_init (&catalog[0].average[Nave]);
    229230    catalog[0].average[Nave].R             = stars[N][0].average.R;
    230231    catalog[0].average[Nave].D             = stars[N][0].average.D;
    231232
    232233    catalog[0].average[Nave].Nmeasure      = NREFSTAR_GROUP;
    233     catalog[0].average[Nave].Nmissing      = 0;
    234     catalog[0].average[Nave].Nextend       = 0;
    235 
    236234    catalog[0].average[Nave].measureOffset = Nmeas;
    237     catalog[0].average[Nave].missingOffset = -1;
    238     catalog[0].average[Nave].extendOffset  = -1;
     235    catalog[0].average[Nave].objID         = objID;
     236    catalog[0].average[Nave].catID         = catID;
     237
     238    if (PSPS_ID) {
     239        catalog[0].average[Nave].extID = CreatePSPSObjectID(catalog[0].average[Nave].R, catalog[0].average[Nave].D);
     240    }
    239241
    240242    if (ACCEPT_MOTION) {
     
    247249      catalog[0].average[Nave].P           = stars[N][0].average.P;
    248250      catalog[0].average[Nave].dP          = stars[N][0].average.dP;
    249     } else {
    250       catalog[0].average[Nave].dR          = 0;
    251       catalog[0].average[Nave].dD          = 0;
    252       catalog[0].average[Nave].uR          = 0;
    253       catalog[0].average[Nave].uD          = 0;
    254       catalog[0].average[Nave].duR         = 0;
    255       catalog[0].average[Nave].duD         = 0;
    256       catalog[0].average[Nave].P           = 0;
    257       catalog[0].average[Nave].dP          = 0;
    258       catalog[0].average[Nave].Xp          = 0;
    259     }
    260 
    261     catalog[0].average[Nave].Xp            = 0;
    262     catalog[0].average[Nave].ChiSqAve      = 0.0;
    263     catalog[0].average[Nave].ChiSqPM       = 0.0;
    264     catalog[0].average[Nave].ChiSqPar      = 0.0;
    265     catalog[0].average[Nave].Tmean         = 0;
    266     catalog[0].average[Nave].Trange        = 0;
    267     catalog[0].average[Nave].Npos          = 0;
    268 
    269     catalog[0].average[Nave].objID         = objID;
    270     catalog[0].average[Nave].catID         = catID;
    271     catalog[0].average[Nave].flags         = 0;
    272     if (PSPS_ID) {
    273         catalog[0].average[Nave].extID = CreatePSPSObjectID(catalog[0].average[Nave].R,
    274                                                             catalog[0].average[Nave].D);
    275     } else {
    276         catalog[0].average[Nave].extID         = 0;
    277     }
    278 
     251    }
    279252
    280253    objID ++;
    281254
    282255    for (j = 0; j < Nsecfilt; j++) {
    283       catalog[0].secfilt[Nave*Nsecfilt+j].M     = NAN;
    284       catalog[0].secfilt[Nave*Nsecfilt+j].dM    = NAN;
    285       catalog[0].secfilt[Nave*Nsecfilt+j].Xm    = NAN_S_SHORT;
    286       catalog[0].secfilt[Nave*Nsecfilt+j].M_20  = NAN_S_SHORT;
    287       catalog[0].secfilt[Nave*Nsecfilt+j].M_80  = NAN_S_SHORT;
    288       catalog[0].secfilt[Nave*Nsecfilt+j].Ncode = 0;
    289       catalog[0].secfilt[Nave*Nsecfilt+j].Nused = 0;
     256      dvo_secfilt_init (&catalog[0].secfilt[Nave*Nsecfilt+j]);
    290257    }
    291258
  • trunk/Ohana/src/addstar/src/mkcmf.c

    r33653 r34260  
    1414void gauss_init (int Nbin);
    1515double rnd_gauss (double mean, double sigma);
     16void writeStars_PS1_V3 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars);
    1617void writeStars_PS1_V2 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars);
    1718void writeStars_PS1_V1 (FTable *ftable, double *X, double *Y, double *M, int Nstars);
     
    276277  if (!strcmp(type, "PS1_V2")) {
    277278    writeStars_PS1_V2 (&ftable, X, Y, M, Flag, Nstars);
     279    found = TRUE;
     280  }
     281  if (!strcmp(type, "PS1_V3")) {
     282    writeStars_PS1_V3 (&ftable, X, Y, M, Flag, Nstars);
    278283    found = TRUE;
    279284  }
     
    572577}
    573578
     579void writeStars_PS1_V3 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars) {
     580
     581  int i;
     582  CMF_PS1_V3 *stars;
     583  float flux, fSN;
     584
     585  // XXX add gaussian-distributed noise based on counts
     586  // this needs to make different output 'stars' entries depending on the desired type
     587  ALLOCATE (stars, CMF_PS1_V3, Nstars);
     588  gauss_init (2048);
     589  for (i = 0; i < Nstars; i++) {
     590    stars[i].detID = i;
     591
     592    flux = pow (10.0, -0.4*M[i]);
     593    fSN = 1.0 / sqrt(flux);
     594
     595    stars[i].X = X[i];
     596    stars[i].Y = Y[i];
     597    stars[i].M = M[i];
     598    stars[i].Map = M[i] - 0.05;
     599
     600    if (ADDNOISE) {
     601      stars[i].X += FX * fSN * rnd_gauss(0.0, 1.0);
     602      stars[i].Y += FY * fSN * rnd_gauss(0.0, 1.0);
     603      stars[i].M += fSN*rnd_gauss(0.0, 1.0);
     604    }
     605
     606    // randomly give poor PSFQF values
     607    if ((BAD_PSFQF_FRAC > 0.0) && (drand48() < BAD_PSFQF_FRAC)) {
     608      stars[i].psfQual   = 0.25;
     609    } else {
     610      stars[i].psfQual   = PSFQUAL;
     611    }
     612   
     613    stars[i].dX = FX * fSN;
     614    stars[i].dY = FY * fSN;
     615    stars[i].dM = fSN;
     616
     617    stars[i].Mpeak     = M[i] + 1.0;
     618    stars[i].sky       = SKY;
     619    stars[i].dSky      = DSKY;
     620    stars[i].psfChisq  = PSFCHI;
     621    stars[i].crNsigma  = CRN;
     622    stars[i].extNsigma = EXTN;
     623    stars[i].fx        = FX;
     624    stars[i].fy        = FY;
     625    stars[i].df        = DF;
     626    stars[i].nFrames   = 1;
     627    stars[i].flags     = Flag[i];
     628
     629    stars[i].kronFlux  = flux * 1.25;
     630    stars[i].kronFluxErr = fSN * flux * 1.25;
     631  }
     632
     633  gfits_table_set_CMF_PS1_V3 (ftable, stars, Nstars);
     634  gfits_modify (ftable->header, "EXTTYPE",   "%s", 1, "PS1_V3");
     635}
     636
  • trunk/Ohana/src/addstar/src/psps_ids.c

    r27527 r34260  
    1919   
    2020uint64_t
     21CreatePSPSStackDetectionID(int sourceID, int imageID, int detID)
     22{
     23  // sourceID : ID of database + table that tracked the image (< 0x100 = 256)
     24  // imageID : external ID of the image which provided the detections (< 0x1000.0000 ~ 2.7e8)
     25  // detID : detection sequence in image (< 0x1000.0000 ~ 2.7e8)
     26
     27  assert (detID    < 0x10000000);
     28  assert (imageID  < 0x10000000);
     29  assert (sourceID < 0x100);
     30 
     31  uint64_t detectid = ((uint64_t)sourceID << 56) + ((uint64_t)imageID << 28) + (uint64_t)detID;
     32  return detectid;
     33}
     34   
     35uint64_t
    2136CreatePSPSObjectID(double ra, double dec)
    2237{
  • trunk/Ohana/src/addstar/src/sky_tessalation.c

    r34088 r34260  
    2323      sky_tessellation_rings (db, level, Nmax);
    2424      return TRUE;
     25    case TAMAS:
     26      sky_tessellation_tamas (db, level, Nmax);
     27      return TRUE;
    2528    default:
    2629      break;
     
    284287    free (ring);
    285288    free (image);
     289  }   
     290  return (TRUE);
     291}
     292
     293// the RINGS tessellation uses the declination zones proposed by Tamas Budavari,
     294// based on code supplied by Tamas 2012.07.23
     295int sky_tessellation_tamas (FITS_DB *db, int level, int Nmax) {
     296
     297  int j, nDEC, Nimage, Nring, Ntotal, Ndigit;
     298  double dec, dDEC;
     299  SkyRectangle *ring;
     300  Image *image;
     301  char format[16];
     302
     303  // The tessellation has one input parameter: the approximate cell size.  Starting with
     304  // the cell size, determine the optimal projection cell height (dDEC) that results in an
     305  // integer number of dec zones between -90 and +90
     306
     307  // in fact, we place a single image on each pole, so the real range of dec is 180.0 - CELLSIZE:
     308
     309  nDEC = (180.0 - CELLSIZE) / CELLSIZE;
     310  dDEC = (180.0 - CELLSIZE) / nDEC;
     311  nDEC += 2;
     312
     313  // how many total projection cells for this realization?  divide sky area by cell area:
     314  // this is used to set the number of digits, so it does not need to be very accurate...
     315  Ntotal = 41254.2 / (dDEC*dDEC);
     316  Ndigit = (int)(log10(Ntotal)) + 1 ;
     317  snprintf (format, 16, "skycell.%%0%dd", Ndigit);
     318
     319  double d2r = M_PI / 180; // is RAD_DEG
     320
     321  // parameter 'a' is the cell size in degrees
     322  double adeg = 3.955;
     323
     324  // half of 'a' in radians and its atan
     325  double halfa = adeg / 2 * d2r;
     326  double halftheta = atan(halfa);
     327
     328  // loop init
     329  dec = 0; // starting Decl. - could change this...
     330 
     331  while (dec < M_PI / 2 - halftheta) {
     332        double dm = dec - halftheta; // eq.5
     333        if (dec == 0) dm = 0; // initial
     334
     335        // dec is modified by the call below
     336        ring = sky_rectangle_tamas (&dec, dm, halfa, halftheta, &Nring, format);
     337        if (!ring) continue;
     338
     339        // subdivide each image (Nx x Ny subcells)
     340        Nimage = NX_SUB*NY_SUB*Nring;
     341        ALLOCATE (image, Image, Nimage);
     342        for (j = 0; j < Nring; j++) {
     343          // convert the SkyRectangles to Images for output
     344          sky_subdivide_image (&image[j*NX_SUB*NY_SUB], &ring[j], NX_SUB, NY_SUB);
     345          // printf("%s %8.2f %8.2f\n", ring[j].name, ring[j].coords.crval1, ring[j].coords.crval2);
     346        }
     347
     348        /* add the new images and save */
     349        dvo_image_addrows (db, image, Nimage);
     350        SetProtect (TRUE);
     351        dvo_image_update (db, VERBOSE);
     352        SetProtect (FALSE);
     353        dvo_image_clear_vtable (db);
     354   
     355        free (ring);
     356        free (image);
    286357  }   
    287358  return (TRUE);
     
    666737}
    667738
     739// define the parameters of a projection centers for this ring
     740// dec : ~ center of ring in Dec
     741// dDEC : approximate height
     742// nring : number of cells generated for this ring
     743// format : guide to generate the filenames (c-type string format)
     744SkyRectangle *sky_rectangle_tamas (double *Dec, double dm, double halfa, double halftheta, int *nring, char *format) {
     745
     746  static int Nname = 0;
     747  int i, j, NX, NY;
     748  SkyRectangle *ring;
     749
     750  double d2r = M_PI / 180; // is RAD_DEG
     751  double dec = *Dec;
     752
     753  int nRA = (int)ceil(M_PI * cos(dm) / halftheta);  // eq.6       
     754  double dRA = 2 * M_PI / nRA; // eq.7
     755  double dp = atan(tan(dec + halftheta) * cos(dRA / 2)); // eq.9
     756
     757  if (dec == 0.0) {
     758    ALLOCATE (ring, SkyRectangle, nRA);
     759  } else {
     760    ALLOCATE (ring, SkyRectangle, 2*nRA);
     761  }
     762
     763  for (i = 0; i < nRA; i++) {
     764    // R.A. can use different phase per ring
     765    double ra = i * dRA; // + phase (watch wraparound)
     766
     767    int npass = (dec == 0.0) ? 1 : 2;
     768    for (j = 0; j < npass; j++) {
     769
     770      int N = j*nRA + i;
     771
     772      memset (&ring[N], 0, sizeof(SkyRectangle));
     773      memset (&ring[N].coords, 0, sizeof(Coords));
     774
     775      ring[N].coords.crval1 = ra / d2r;
     776      ring[N].coords.crval2 = (j == 0) ? dec / d2r : -dec / d2r;
     777
     778      printf(" \t %d   %25.20f   %25.20f\n", i, ring[N].coords.crval2, ring[N].coords.crval1);
     779
     780      ring[N].coords.pc1_1 = +1.0 * X_PARITY;
     781      ring[N].coords.pc1_2 = +0.0;
     782      ring[N].coords.pc2_1 = -0.0;
     783      ring[N].coords.pc2_2 = +1.0;
     784 
     785      // range values are in projected degrees
     786      NX = cos(dec - halftheta) * dRA   * 3600.0 / SCALE / d2r;
     787      NY =    2 * halftheta * 3600.0 / SCALE / d2r;
     788
     789      // crpix1,crpix2 is the projection center
     790      ring[N].coords.crpix1 = 0.5*NX;
     791      ring[N].coords.crpix2 = 0.5*NY;
     792
     793      ring[N].coords.cdelt1 = SCALE / 3600.0;
     794      ring[N].coords.cdelt2 = SCALE / 3600.0;
     795
     796      strcpy (ring[N].coords.ctype, "DEC--TAN");
     797
     798      ring[N].NX = NX*(1.0 + PADDING);
     799      ring[N].NY = NY*(1.0 + PADDING);
     800      ring[N].photcode = 1; // this needs to be set more sensibly
     801
     802      snprintf (ring[N].name, DVO_IMAGE_NAME_LEN, format, Nname);
     803      Nname++;
     804    }
     805  }
     806
     807  // advance to next ring
     808  *Dec = halftheta + dp;
     809
     810  *nring = (dec == 0.0) ? nRA : 2*nRA;
     811  return ring;
     812}
     813
    668814// an allocated image set is supplied, we fill in the values
    669815int sky_subdivide_image (Image *output, SkyRectangle *input, int Nx, int Ny) {
     
    682828  }
    683829
    684   Ndigit = (int)(log10(Nx*Ny)) + 1 ;
    685   snprintf (format, 24, "%s.%%0%dd", input[0].name, Ndigit);
     830  if (Nx * Ny > 1) {
     831    Ndigit = (int)(log10(Nx*Ny)) + 1 ;
     832    snprintf (format, 24, "%s.%%0%dd", input[0].name, Ndigit);
     833  } else {
     834    snprintf (format, 24, "%s", input[0].name);
     835  }
    686836
    687837  // if requested extend, the skycell boundaries so that skycells overlap
     
    696846      memcpy (&output[N].coords, &input[0].coords, sizeof(Coords));
    697847
    698       snprintf (output[N].name, DVO_IMAGE_NAME_LEN, format, N);
     848      if (Nx + Ny > 1) {
     849        snprintf (output[N].name, DVO_IMAGE_NAME_LEN, format, N);
     850      } else {
     851        snprintf (output[N].name, DVO_IMAGE_NAME_LEN, "%s", format);
     852      }
     853
    699854      output[N].NX = NX + 2 * pad_x;
    700855      output[N].NY = NY + 2 * pad_y;
  • trunk/Ohana/src/addstar/test/simple.dvo

    r33653 r34260  
    2121  test.fields PS1_V2    PS1_V3
    2222  test.fields PS1_V3    PS1_V3
     23
     24  test.fields PS1_DEV_0 PS1_V4
     25  test.fields PS1_DEV_1 PS1_V4
     26  test.fields PS1_V1    PS1_V4
     27  test.fields PS1_V2    PS1_V4
     28  test.fields PS1_V3    PS1_V4
    2329end 
    2430
     
    8389    sort id1 v1
    8490    sort id2 v2
     91
     92    # some fields require arithmetic manipulations
     93    if ("$name:0" == "KRON_FLUX")
     94     set v1 = -2.5*log(v1)
     95    end
     96    if ("$name:0" == "KRON_FLUX_ERR")
     97     set v1 = KRON_FLUX_ERR / KRON_FLUX
     98    end
     99
    85100    set d = v1 - v2
    86101    vstat -q d
     
    88103    #echo tapOK fabs($MEAN)  < 0.001 "$name:0 vs $name:2 (MEAN)"
    89104    #echo tapOK fabs($SIGMA) < 0.001 "$name:0 vs $name:2 (SIGMA)"
     105
     106    # THETA is stored to only (360/65536) deg accuracy
     107    if ("$name:0" == "PSF_THETA")
     108      echo $MEAN
     109      tapOK {abs($MEAN)  < 0.006} "$name:0 vs $name:2 (MEAN)"
     110      tapOK {abs($SIGMA) < 0.001} "$name:0 vs $name:2 (SIGMA)"
     111      continue
     112    end
    90113
    91114    tapOK {abs($MEAN)  < 0.001} "$name:0 vs $name:2 (MEAN)"
     
    111134  output stdout
    112135end
     136
     137# the following lists define fields in the cmf files which can be compared to their equivalents in DVO
     138# the left column is the cmf field name, the right column is the dvo field name
    113139
    114140# list of cmf fields to test matched to mextract fields
     
    122148  PSF_INST_MAG      : mag:inst
    123149  PSF_INST_MAG_SIG  : mag:err
    124   PEAK_FLUX_AS_MAG  : SKIP
     150  PEAK_FLUX_AS_MAG  : SKIP # not ingested into DVO
    125151  SKY               : sky
    126152  SKY_SIG           : sky_err
     
    130156  PSF_THETA         : THETA
    131157  PSF_QF            : PSF_QF
    132   N_FRAMES          : SKIP
     158  N_FRAMES          : SKIP # not ingested into DVO
    133159end
    134160
     
    143169  PSF_INST_MAG      : mag:inst
    144170  PSF_INST_MAG_SIG  : mag:err
    145   PEAK_FLUX_AS_MAG  : SKIP
     171  PEAK_FLUX_AS_MAG  : SKIP # not ingested into DVO
    146172  SKY               : sky
    147173  SKY_SIG           : sky_err
     
    153179  PSF_THETA         : THETA
    154180  PSF_QF            : PSF_QF
    155   N_FRAMES          : SKIP
     181  N_FRAMES          : SKIP # not ingested into DVO
    156182  FLAGS             : phot_flags
    157183end
     
    236262  X_PSF_SIG         : xccd:err # FAIL
    237263  Y_PSF_SIG         : yccd:err # FAIL
    238   RA_PSF            : SKIP # astrometry is not calibrated in the cmf
    239   DEC_PSF           : SKIP # astrometry is not calibrated in the cmf
    240264  POSANGLE          : SKIP # astrometry is not calibrated in the cmf
    241265  PLTSCALE          : SKIP # astrometry is not calibrated in the cmf
    242266  PSF_INST_MAG      : mag:inst 
    243267  PSF_INST_MAG_SIG  : mag:err   
    244   AP_MAG_STANDARD   : mag:ap # FAIL
    245   AP_MAG_RADIUS     : SKIP # no accessor
    246   PEAK_FLUX_AS_MAG  : SKIP # no accessor
     268  PSF_INST_FLUX     : SKIP # not ingested into DVO
     269  PSF_INST_FLUX_SIG : SKIP # not ingested into DVO
     270  AP_MAG_STANDARD   : mag:aperinst # FAIL
     271  AP_MAG_RAW        : SKIP # not ingested into DVO
     272  AP_MAG_RADIUS     : SKIP # not ingested into DVO
    247273  CAL_PSF_MAG       : SKIP # photometry is not calibrated in the cmf
    248274  CAL_PSF_MAG_SIG   : SKIP # photometry is not calibrated in the cmf
     275  RA_PSF            : SKIP # astrometry is not calibrated in the cmf
     276  DEC_PSF           : SKIP # astrometry is not calibrated in the cmf
     277  PEAK_FLUX_AS_MAG  : SKIP # not ingested into DVO
    249278  SKY               : sky       
    250279  SKY_SIG           : sky_err   
     
    256285  PSF_THETA         : THETA # FAIL
    257286  PSF_QF            : PSF_QF   
    258   PSF_NDOF          : SKIP # no accessor
    259   PSF_NPIX          : SKIP # no accessor
    260   MOMENTS_XX        : SKIP # no accessor
    261   MOMENTS_XY        : SKIP # no accessor
    262   MOMENTS_YY        : SKIP # no accessor
     287  PSF_QF_PERFECT    : SKIP # not ingested into DVO
     288  PSF_NDOF          : PSF_NDOF
     289  PSF_NPIX          : PSF_NPIX
     290  MOMENTS_XX        : MXX
     291  MOMENTS_XY        : MXY
     292  MOMENTS_YY        : MYY
     293  MOMENTS_M3C       : SKIP # not ingested into DVO
     294  MOMENTS_M3S       : SKIP # not ingested into DVO
     295  MOMENTS_M4C       : SKIP # not ingested into DVO
     296  MOMENTS_M4S       : SKIP # not ingested into DVO
     297  MOMENTS_R1        : SKIP # not ingested into DVO
     298  MOMENTS_RH        : SKIP # not ingested into DVO
     299  KRON_FLUX         : mag:kroninst
     300  KRON_FLUX_ERR     : mag:kronerr
     301  KRON_FLUX_INNER   : SKIP # not ingested into DVO
     302  KRON_FLUX_OUTER   : SKIP # not ingested into DVO
    263303  FLAGS             : phot_flags
    264   N_FRAMES          : SKIP # no accessor       
    265 end
     304  N_FRAMES          : SKIP # not ingested into DVO
     305end
Note: See TracChangeset for help on using the changeset viewer.