IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 23941


Ignore:
Timestamp:
Apr 21, 2009, 1:34:06 AM (17 years ago)
Author:
Sebastian Jester
Message:

Insert comments with proposals/questions; prepare to make changes to code to filter out objects that are not OK_RUN

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/sj_branches/sj_SDSSaddstar_branch_r23928_20090420/Ohana/src/addstar/src/ReadStarsSDSS.c

    r21508 r23941  
    2929// XXX NOTE : as of 2008.02.27, the zero point is still carried internally in millimags
    3030
    31 // given a file with the pointer at the start of the table block and the 
     31// given a file with the pointer at the start of the table block and the
    3232// corresponding image header, load the stars from the table
    3333Stars *ReadStarsSDSS (FILE *f, char *name, Header *header, Header *in_theader, Image *images, int *nimages, unsigned int *nstars) {
     
    4444  int photcode[5];
    4545  int Nrow, Ncol; // used in the GET_COLUMN_1,5 macros above
    46  
     46
    4747  if (in_theader == NULL) {
    4848    table.header = &theader;
     
    5151    table.header = in_theader;
    5252    Nskip = in_theader[0].size;
    53     fseek (f, Nskip, SEEK_CUR); 
     53    fseek (f, Nskip, SEEK_CUR);
    5454  }
    5555
     
    6666  strcpy (filtname[4], "z");
    6767
     68  // XXXYYY Why these? and not the 'ref' version? see 'Why these?' below, too.
    6869  NAMED_PHOTCODE_AND_ZP (photcode[0], zeropt[0], "U_SDSS");
    6970  NAMED_PHOTCODE_AND_ZP (photcode[1], zeropt[1], "G_SDSS");
     
    9293  tzero[4] = ohana_mjd_to_sec (mjd[4]);
    9394
     95  // XXXYYY EDR paper says "psfWidth is is the effective width also determined at the center of each frame. It
     96  // is a good generic number to quote for the seeing on each frame." It is stored in the tsField file as
     97  // psf_width.
    9498  gfits_scan (table.header, "SEEING_U", "%f", 1, &seeing[0]);
    9599  gfits_scan (table.header, "SEEING_G", "%f", 1, &seeing[1]);
     
    98102  gfits_scan (table.header, "SEEING_Z", "%f", 1, &seeing[4]);
    99103
     104  // XXXYYY What's the meaning of photErr? Where is it used again?
    100105  gfits_scan (table.header, "PSFERR_U", "%f", 1, &photErr[0]);
    101106  gfits_scan (table.header, "PSFERR_G", "%f", 1, &photErr[1]);
     
    106111  gfits_scan (header, "CAMCOL", "%d", 1, &camcol); // value in header is usec / unbinned row
    107112
     113  // XXXYYY which number is read here? In other words, where is the zeropoint set that is returned by
     114  // GetZeroPoint()? I can't find the corresponding SetZeroPoint()
    108115  ZeroPt = GetZeroPoint();
    109116
     
    132139  GET_COLUMN_5 (psfCountsErr, float);
    133140
     141  GET_COLUMN_1 (status, int);
     142
    134143  // the value of stars[].M is supposed to be the instrumental magnitude offset by the
    135144  // default zero point 25.0 (-2.5*log_10(counts/sec) + ZeroPt).  The magnitude reported
     
    141150      N = NFILTER*i + j;
    142151      InitStar (&stars[N]);
    143      
     152
    144153      // any values not explicitly set are left at 0.0
    145154      stars[N].average.R         = ra[i] + dCOS(dec[i]) * offsetRa[N] / 3600.0;
     
    152161      stars[N].measure.dXccd     = ShortPixels(colcErr[N]);
    153162      stars[N].measure.dYccd     = ShortPixels(rowcErr[N]);
     163      // XXXYYY again the question what ZeroPt-zeropt[j] evalutes to - according to what it says above, this
     164      // should take .M back into instrumental magnitudes, but does it? In particular, the *actual* zeropoint
     165      // that was used for the computation of this calibrated magnitude is in the tsField file; in general, it
     166      // will differ by some small amount (less than 1% I'm guessing) from the 'default' zero points; for
     167      // non-photometric stripe82 data, the actual zero point will differ A LOT from the default one.
    154168      stars[N].measure.M         = psfCounts[N] + ZeroPt - zeropt[j];
    155169      stars[N].measure.dM        = psfCountsErr[N];
    156170      stars[N].measure.Map       = fiberCounts[N] + ZeroPt - zeropt[j];
    157       stars[N].measure.Sky       = sky[N]; // adjust this to counts?
     171      stars[N].measure.Sky       = sky[N]; // adjust this to counts? XXXYYY It's given in maggies; to get back to counts,
    158172      stars[N].measure.dSky      = skyErr[N];
    159       stars[N].measure.FWx       = ShortPixels(seeing[j]); // reported in arcsec?
     173      stars[N].measure.FWx       = ShortPixels(seeing[j]); // reported in arcsec? XXXYYY: Yes!
    160174      stars[N].measure.FWy       = ShortPixels(seeing[j]);
    161175      stars[N].measure.psfChisq  = prob_psf[N]; // XXX not really the correct value...
    162176      stars[N].measure.detID     = N;
    163177      stars[N].measure.t         = tzero[j] + clockRate*rowc[N]; // time since row 0
    164       stars[N].measure.dt        = 53.907456; // is this 2048*clockRate ?
     178      stars[N].measure.dt        = 53.907456; // is this 2048*clockRate ? XXXYYY Nearly - 4 milliseconds more. shouldn't this be a #define?
    165179
    166180      SetSDSSFlags (&stars[N], flags[N], flags2[N]);
     
    178192      altaz (&alt, &az, sidtime - stars[N].average.R, stars[N].average.D, Latitude);
    179193
     194      // XXXYYY airmass is in tsField file
    180195      stars[N].measure.airmass   = 1.0 / dCOS(90.0 - alt);
    181196      stars[N].measure.az        = az;
    182197      stars[N].measure.photcode  = photcode[j];
    183198    }
    184   }   
     199  }
    185200
    186201  for (i = 0; i < NFILTER; i++) {
    187202
    188203    N = i + *nimages;
    189    
     204
    190205    // XXX for now, we define a totally fake coordinate system centered on the first listed star
     206
     207    // XXXYYY what's happening here? *Exact* astrometry for SDSS can only be retrieved using the
     208    // colour-dependent distortion coefficients from the tsField file. See
     209    // www.sdss.org/dr7/products/general/astrometry.html
     210    //
     211    // WCS headers are of course present in the fpC image,
     212    // but then, that would be downloading a lot more data.
    191213    strcpy (images[N].coords.ctype, "RA---TAN");
    192    
     214
    193215    images[N].coords.crval1 = stars[0].average.R;
    194216    images[N].coords.crval2 = stars[0].average.D;
    195    
     217
    196218    images[N].coords.crpix1 = stars[0].measure.Xccd;
    197219    images[N].coords.crpix2 = stars[0].measure.Yccd;
     
    207229
    208230    images[N].NX = 2048;
     231    // XXXYYY NZ should be 1361, since the last 128 rows are the duplicate of the neighbouring frame, and only
     232    // objects in 64 <= row < 1361+64 (or maybe the <= is at the upper limit?) are counted as belonging to
     233    // "this" field.  I.e. we also need to read the status flag and use only primary and secondary objects
     234    // (i.e. those that have OK_RUN set in status); otherwise, the SAME PHOTONS will appear more than once in
     235    // the dvo database.
     236
     237    // From http://www.sdss.org/DR7/dm/flatFiles/tsObj.html#status
     238
     239    // AR_OBJECT_STATUS_OK_RUN = 0x10, /* Located within the primary range */
     240                                        /* of rows for this field.  This is  */
     241                                        /* usable object.  This flag is set  */
     242                                        /* by "setObjectStatus".             */
     243    // OK_RUN is the union of the following two classes:
     244    // AR_OBJECT_STATUS_SECONDARY  = 0x1000, /* This is a secondary survey object.*/
     245    // AR_OBJECT_STATUS_PRIMARY    = 0x2000, /* This is a primary survey object.  */    // 'primary' means 'inside the area that defines a unique coverage of the sky', i.e. selecting primary
     246    // objects rejects repeat observations by design (the SDSS philosophy: everything gets observed
     247    // once). Adding 'secondary' objects includes repeat observations of the same object from a different
     248    // run. This is desired in the case of dvo. Hence we need exactly those objects with OK_RUN set, not more,
     249    // not less.
     250
     251    // And once more, NY  = 1361, which led into the whole 'status' discussion.
    209252    images[N].NY = 1490;
    210253
    211254    images[N].tzero = tzero[i];
    212255    images[N].cerror = 0.0;
    213  
     256
    214257    // set photcodes for the 5 images (SDSS_U,G,R,I,Z)
     258    // XXXYYY Why these? and not the 'u_sdss etc. 'ref' photcodes (see above where it says U_SDSS etc.)
    215259    images[N].photcode = photcode[i];
    216260
    217261    // calculate this from : C_OBS, TRACKING, and NY
    218262    images[N].exptime = 2048*clockRate;
    219  
     263
    220264    images[N].apmifit = 0.0;
    221265    images[N].dapmifit = 0.0;
    222     images[N].detection_limit = 0.0; 
     266    images[N].detection_limit = 0.0;
    223267    images[N].saturation_limit = 0.0;
    224268    images[N].fwhm_x = seeing[i];
     
    249293
    250294    images[N].nstar = Nstars;
    251  
     295
    252296    images[N].imageID = N;
    253297    images[N].externID = 0;
     
    285329int SetSDSSFlags (Stars *star, unsigned int flags1, unsigned int flags2) {
    286330
    287   // XXX this is wrong, need to roll left to set the correct bit 
     331  // XXX this is wrong, need to roll left to set the correct bit
    288332  if (flags1 & 0x00000002) star[0].measure.photFlags |= 0x0001; // BRIGHT            - 1  1
    289333  if (flags1 & 0x00000004) star[0].measure.photFlags |= 0x0002; // EDGE              - 1  2
Note: See TracChangeset for help on using the changeset viewer.