IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 5585


Ignore:
Timestamp:
Nov 23, 2005, 10:47:32 AM (20 years ago)
Author:
eugene
Message:

fixed up calibration process

Files:
9 edited

Legend:

Unmodified
Added
Removed
  • branches/addstar-1-3-0/Ohana/src/addstar/include/addstar.h

    r5347 r5585  
    6868/* these globals are used separately by both client and server (KEEP) */
    6969double ZeroPt;  // double check for consistency
     70double CAL_INSTMAG_MAX;
     71double CAL_INSTMAG_MIN;
    7072int    VERBOSE;
    7173
     
    9395void       InitCalibration        PROTO(());
    9496GSCRegion *LoadRegions            PROTO((int *nregions));
    95 void       SaveCalibration        PROTO((float M, float dM, float Mr, float dMr, float Mc, float A, int N));
     97void       SaveCalibration        PROTO((float Mo, float dMo, float Mr, float dMr, float Mi, int N));
    9698void       SetProtect             PROTO((int mode));
    9799int        SetSignals             PROTO(());
  • branches/addstar-1-3-0/Ohana/src/addstar/src/ConfigInit.c

    r5347 r5585  
    4343  ScanConfig (config, "OBSERVATORY-LATITUDE",   "%lf", 0, &Latitude);
    4444  ScanConfig (config, "SUBPIX_DATAFILE",        "%s",  0, SubpixDatafile);
     45
     46  /* instrumental magnitude range for calibration mode */
     47  if (!ScanConfig (config, "CAL_INSTMAG_MAX",   "%lf", 0, &CAL_INSTMAG_MAX)) CAL_INSTMAG_MAX =  -9.0;
     48  if (!ScanConfig (config, "CAL_INSTMAG_MIN",   "%lf", 0, &CAL_INSTMAG_MIN)) CAL_INSTMAG_MIN = -13.0;
    4549
    4650  /* location of needed data sources */
  • branches/addstar-1-3-0/Ohana/src/addstar/src/args.c

    r5328 r5585  
    144144    options.calibrate = TRUE;
    145145    remove_argument (N, &argc, argv);
     146    InitCalibration (FALSE);
     147  }
     148  if ((N = get_argument (argc, argv, "-excal"))) {
     149    options.calibrate = TRUE;
     150    remove_argument (N, &argc, argv);
     151    InitCalibration (FALSE);
     152  }
     153  if ((N = get_argument (argc, argv, "-incal"))) {
     154    options.calibrate = TRUE;
     155    remove_argument (N, &argc, argv);
     156    InitCalibration (TRUE);
    146157  }
    147158
  • branches/addstar-1-3-0/Ohana/src/addstar/src/calibrate.c

    r5328 r5585  
    11# include "addstar.h"
    22
     3static int InternalCal;
    34static int Ncal, NCAL, *Nstar;
    4 static float *Mobs, *dMobs, *Mref, *dMref, *Cref, *Airm;
     5static float *Mobs, *dMobs, *Mref, *dMref, *Minst;
    56
    6 void InitCalibration () {
     7void InitCalibration (int mode) {
    78
    8     fprintf (stderr, "calibrating the image...\n");
    9     Ncal = 0;
    10     NCAL = 1000;
    11     ALLOCATE (Mobs,  float,  NCAL);
    12     ALLOCATE (dMobs, float, NCAL);
    13     ALLOCATE (Mref,  float,  NCAL);
    14     ALLOCATE (dMref,  float,  NCAL);
    15     ALLOCATE (Cref,  float,  NCAL);
    16     ALLOCATE (Airm,  float,  NCAL);
    17     ALLOCATE (Nstar, int, NCAL);
     9  InternalCal = mode;
     10
     11  fprintf (stderr, "calibrating the image...\n");
     12  Ncal = 0;
     13  NCAL = 1000;
     14  ALLOCATE (Mobs,  float,  NCAL);
     15  ALLOCATE (dMobs, float, NCAL);
     16  ALLOCATE (Mref,  float,  NCAL);
     17  ALLOCATE (dMref, float,  NCAL);
     18  ALLOCATE (Minst, float,  NCAL);
     19  ALLOCATE (Nstar, int, NCAL);
    1820}
    1921 
    20 void SaveCalibration (float M, float dM, float Mr, float dMr, float Mc, float A, int N) {
     22void SaveCalibration (float Mo, float dMo, float Mr, float dMr, float Mi, int N) {
    2123
    22   Mobs[Ncal] = M;
    23   dMobs[Ncal] = dM;
    24   Mref[Ncal] = Mr;
     24  Mobs[Ncal]  = Mo;
     25  dMobs[Ncal] = dMo;
     26  Mref[Ncal]  = Mr;
    2527  dMref[Ncal] = dMr;
    26   Cref[Ncal] = Mc;
    27   Airm[Ncal] = A;
     28  Minst[Ncal] = Mi;
    2829  Nstar[Ncal] = N;
    2930  Ncal ++;
     
    3132  if (Ncal == NCAL) {
    3233    NCAL += 1000;
     34    REALLOCATE (Mobs,  float, NCAL);
    3335    REALLOCATE (dMobs, float, NCAL);
    34     REALLOCATE (Mobs,  float, NCAL);
    3536    REALLOCATE (Mref,  float, NCAL);
    36     REALLOCATE (dMref,  float, NCAL);
    37     REALLOCATE (Cref,  float, NCAL);
    38     REALLOCATE (Airm,  float, NCAL);
     37    REALLOCATE (dMref, float, NCAL);
     38    REALLOCATE (Minst, float, NCAL);
    3939    REALLOCATE (Nstar, int,   NCAL);
    4040  }
     
    4444void AddToCalibration (Average *average, Measure *measure, Measure *new, int *next, int Nstar) {
    4545
    46   int i, m, found0, found1, found2;
     46  int i, j, m, found0, found1, found2;
    4747  float CalM0, CalM1, CalM2, dCalM;
     48  float Mcal, color, factor;
    4849  short CalC0, CalC1, CalC2;
    49   PhotCode *code;
     50
     51  PhotCode *mycode;  // photcode of this measurement
     52  PhotCode *incode;  // mycode.equiv (internal reference)
     53  PhotCode *excode;  // incode.equiv (external reference)
    5054
    5155  found0 = found1 = found2 = FALSE;
    5256  CalM0 = CalM1 = CalM2 = dCalM = NO_MAG;
    5357
    54   code  = GetPhotcodebyCode (new[0].source);
    55   CalC0 = code[0].equiv;
    56   CalC1 = code[0].c1;
    57   CalC2 = code[0].c2;
     58  // we have two options here:
     59  //  - calibrate to internal system (Mcal)
     60  //  - calibrate to external system (Mref)
     61
     62  mycode = GetPhotcodebyCode (new[0].source);
     63  incode = GetPhotcodebyCode (mycode[0].equiv);
     64  excode = GetPhotcodebyCode (incode[0].equiv);
     65
     66  if (InternalCal) {
     67    CalC0 = incode[0].code;
     68  } else {
     69    CalC0 = excode[0].code;
     70  }
     71  CalC1 = mycode[0].c1;
     72  CalC2 = mycode[0].c2;
    5873
    5974  m = average[0].offset;
     
    7388    }
    7489    if (found0 && found1 && found2) {
    75       SaveCalibration (new[0].M_PS, new[0].dM_PS, CalM0, dCalM, CalM1-CalM2, new[0].airmass_PS, Nstar);
     90      Mcal   = new[0].M_PS + 0.001*mycode[0].C + mycode[0].K*(new[0].airmass_PS - 1.0) - ZeroPt;
     91      color  = CalM1 - CalM2 - 0.001*mycode[0].dX;
     92      factor = color;
     93      for (j = 0; j < mycode[0].Nc; j++) {
     94        Mcal += mycode[0].X[j]*factor;
     95        factor *= color;
     96      }
     97      if (!InternalCal) {
     98        Mcal += 0.001*incode[0].C;
     99      }
     100      // if we want to apply a Mcal -> Mref color correction, we need the additional color term
     101      SaveCalibration (Mcal, new[0].dM_PS, CalM0, dCalM, new[0].M_PS - ZeroPt - new[0].dt_PS, Nstar);
    76102      return;
    77103    }
     
    84110
    85111  int i, MaxN, *Nlist, Nkeep;
    86   float Mint;
    87   float N, M1, M2, Klam, Clam, Xlam, Mabs, *Dmag, *dDmag;
     112  float N, M1, M2, *Dmag, *dDmag;
    88113  float dMo, dMr, Mw, Dmed, W1, W2, NSigma;
    89   PhotCode *code;
    90 
    91   code = GetPhotcodebyCode (image[0].source);
    92114
    93115  /* reject multiple matched-stars */
     
    114136  ALLOCATE (dDmag, float, Ncal);
    115137  Nkeep = 0;
    116   Clam = code[0].C*0.001; /* photcode.C still in millimags */
    117   Klam = code[0].K;
    118   Xlam = code[0].X[0];
    119138  for (i = 0; i < Ncal; i++) {
    120139    /* if this entry has too many (or two few?) matches, skip it */
    121140    if (Nlist[Nstar[i]] != 1) continue;
    122     Mint = Clam + Mobs[i] + Xlam*Cref[i] + Klam*(Airm[i] - 1.000);
    123     Mabs = Mint - ZeroPt;
    124     /* note: subpix correction is applied in gstars */
    125141
    126     /* skip stars brighter than 8.0 */
    127     /* XXX EAM : perhaps make this more flexible?? */
    128     if (Mabs > 9.0) continue;
    129     if (Mabs < 5.0) continue;
     142    /* clip by instrumental magnitude */
     143    if (Minst[i] > CAL_INSTMAG_MAX) continue;
     144    if (Minst[i] < CAL_INSTMAG_MIN) continue;
    130145   
    131146    /* XXX EAM: note the artificial 0.005 dmag here */
     
    133148    dMo = MAX (0.005, dMobs[i]);
    134149
    135     Dmag[Nkeep] = (Mabs - Mref[i]);
     150    Dmag[Nkeep] = (Mobs[i] - Mref[i]);
    136151    dDmag[Nkeep] = (dMr*dMr + dMo*dMo);
    137152    Nkeep ++;
  • trunk/Ohana/src/addstar/include/addstar.h

    r5445 r5585  
    7171/* these globals are used separately by both client and server (KEEP) */
    7272double ZeroPt;  // double check for consistency
     73double CAL_INSTMAG_MAX;
     74double CAL_INSTMAG_MIN;
    7375int    VERBOSE;
    7476
     
    9193FILE      *GetDB                  PROTO((int *state));
    9294void       InitCalibration        PROTO(());
    93 void       SaveCalibration        PROTO((float M, float dM, float Mr, float dMr, float Mc, float A, int N));
     95void       SaveCalibration        PROTO((float Mo, float dMo, float Mr, float dMr, float Mi, int N));
    9496void       SetProtect             PROTO((int mode));
    9597int        SetSignals             PROTO(());
  • trunk/Ohana/src/addstar/src/ConfigInit.c

    r5445 r5585  
    4343  ScanConfig (config, "OBSERVATORY-LATITUDE",   "%lf", 0, &Latitude);
    4444  ScanConfig (config, "SUBPIX_DATAFILE",        "%s",  0, SubpixDatafile);
     45
     46  /* instrumental magnitude range for calibration mode */
     47  if (!ScanConfig (config, "CAL_INSTMAG_MAX",   "%lf", 0, &CAL_INSTMAG_MAX)) CAL_INSTMAG_MAX =  -9.0;
     48  if (!ScanConfig (config, "CAL_INSTMAG_MIN",   "%lf", 0, &CAL_INSTMAG_MIN)) CAL_INSTMAG_MIN = -13.0;
    4549
    4650  /* location of needed data sources */
  • trunk/Ohana/src/addstar/src/addstar.c

    r5448 r5585  
    6060    SkyList *tmp;
    6161    tmp = SkyListExistingSubset (skylist, CATDIR);
    62     SkyListFree (skylist);
     62    SkyListFree (skylist, FALSE);
    6363    skylist = tmp;
    6464  }
  • trunk/Ohana/src/addstar/src/args.c

    r5443 r5585  
    144144    options.calibrate = TRUE;
    145145    remove_argument (N, &argc, argv);
     146    InitCalibration (FALSE);
     147  }
     148  if ((N = get_argument (argc, argv, "-excal"))) {
     149    options.calibrate = TRUE;
     150    remove_argument (N, &argc, argv);
     151    InitCalibration (FALSE);
     152  }
     153  if ((N = get_argument (argc, argv, "-incal"))) {
     154    options.calibrate = TRUE;
     155    remove_argument (N, &argc, argv);
     156    InitCalibration (TRUE);
    146157  }
    147158
  • trunk/Ohana/src/addstar/src/calibrate.c

    r5328 r5585  
    11# include "addstar.h"
    22
     3static int InternalCal;
    34static int Ncal, NCAL, *Nstar;
    4 static float *Mobs, *dMobs, *Mref, *dMref, *Cref, *Airm;
     5static float *Mobs, *dMobs, *Mref, *dMref, *Minst;
    56
    6 void InitCalibration () {
     7void InitCalibration (int mode) {
    78
    8     fprintf (stderr, "calibrating the image...\n");
    9     Ncal = 0;
    10     NCAL = 1000;
    11     ALLOCATE (Mobs,  float,  NCAL);
    12     ALLOCATE (dMobs, float, NCAL);
    13     ALLOCATE (Mref,  float,  NCAL);
    14     ALLOCATE (dMref,  float,  NCAL);
    15     ALLOCATE (Cref,  float,  NCAL);
    16     ALLOCATE (Airm,  float,  NCAL);
    17     ALLOCATE (Nstar, int, NCAL);
     9  InternalCal = mode;
     10
     11  fprintf (stderr, "calibrating the image...\n");
     12  Ncal = 0;
     13  NCAL = 1000;
     14  ALLOCATE (Mobs,  float,  NCAL);
     15  ALLOCATE (dMobs, float, NCAL);
     16  ALLOCATE (Mref,  float,  NCAL);
     17  ALLOCATE (dMref, float,  NCAL);
     18  ALLOCATE (Minst, float,  NCAL);
     19  ALLOCATE (Nstar, int, NCAL);
    1820}
    1921 
    20 void SaveCalibration (float M, float dM, float Mr, float dMr, float Mc, float A, int N) {
     22void SaveCalibration (float Mo, float dMo, float Mr, float dMr, float Mi, int N) {
    2123
    22   Mobs[Ncal] = M;
    23   dMobs[Ncal] = dM;
    24   Mref[Ncal] = Mr;
     24  Mobs[Ncal]  = Mo;
     25  dMobs[Ncal] = dMo;
     26  Mref[Ncal]  = Mr;
    2527  dMref[Ncal] = dMr;
    26   Cref[Ncal] = Mc;
    27   Airm[Ncal] = A;
     28  Minst[Ncal] = Mi;
    2829  Nstar[Ncal] = N;
    2930  Ncal ++;
     
    3132  if (Ncal == NCAL) {
    3233    NCAL += 1000;
     34    REALLOCATE (Mobs,  float, NCAL);
    3335    REALLOCATE (dMobs, float, NCAL);
    34     REALLOCATE (Mobs,  float, NCAL);
    3536    REALLOCATE (Mref,  float, NCAL);
    36     REALLOCATE (dMref,  float, NCAL);
    37     REALLOCATE (Cref,  float, NCAL);
    38     REALLOCATE (Airm,  float, NCAL);
     37    REALLOCATE (dMref, float, NCAL);
     38    REALLOCATE (Minst, float, NCAL);
    3939    REALLOCATE (Nstar, int,   NCAL);
    4040  }
     
    4444void AddToCalibration (Average *average, Measure *measure, Measure *new, int *next, int Nstar) {
    4545
    46   int i, m, found0, found1, found2;
     46  int i, j, m, found0, found1, found2;
    4747  float CalM0, CalM1, CalM2, dCalM;
     48  float Mcal, color, factor;
    4849  short CalC0, CalC1, CalC2;
    49   PhotCode *code;
     50
     51  PhotCode *mycode;  // photcode of this measurement
     52  PhotCode *incode;  // mycode.equiv (internal reference)
     53  PhotCode *excode;  // incode.equiv (external reference)
    5054
    5155  found0 = found1 = found2 = FALSE;
    5256  CalM0 = CalM1 = CalM2 = dCalM = NO_MAG;
    5357
    54   code  = GetPhotcodebyCode (new[0].source);
    55   CalC0 = code[0].equiv;
    56   CalC1 = code[0].c1;
    57   CalC2 = code[0].c2;
     58  // we have two options here:
     59  //  - calibrate to internal system (Mcal)
     60  //  - calibrate to external system (Mref)
     61
     62  mycode = GetPhotcodebyCode (new[0].source);
     63  incode = GetPhotcodebyCode (mycode[0].equiv);
     64  excode = GetPhotcodebyCode (incode[0].equiv);
     65
     66  if (InternalCal) {
     67    CalC0 = incode[0].code;
     68  } else {
     69    CalC0 = excode[0].code;
     70  }
     71  CalC1 = mycode[0].c1;
     72  CalC2 = mycode[0].c2;
    5873
    5974  m = average[0].offset;
     
    7388    }
    7489    if (found0 && found1 && found2) {
    75       SaveCalibration (new[0].M_PS, new[0].dM_PS, CalM0, dCalM, CalM1-CalM2, new[0].airmass_PS, Nstar);
     90      Mcal   = new[0].M_PS + 0.001*mycode[0].C + mycode[0].K*(new[0].airmass_PS - 1.0) - ZeroPt;
     91      color  = CalM1 - CalM2 - 0.001*mycode[0].dX;
     92      factor = color;
     93      for (j = 0; j < mycode[0].Nc; j++) {
     94        Mcal += mycode[0].X[j]*factor;
     95        factor *= color;
     96      }
     97      if (!InternalCal) {
     98        Mcal += 0.001*incode[0].C;
     99      }
     100      // if we want to apply a Mcal -> Mref color correction, we need the additional color term
     101      SaveCalibration (Mcal, new[0].dM_PS, CalM0, dCalM, new[0].M_PS - ZeroPt - new[0].dt_PS, Nstar);
    76102      return;
    77103    }
     
    84110
    85111  int i, MaxN, *Nlist, Nkeep;
    86   float Mint;
    87   float N, M1, M2, Klam, Clam, Xlam, Mabs, *Dmag, *dDmag;
     112  float N, M1, M2, *Dmag, *dDmag;
    88113  float dMo, dMr, Mw, Dmed, W1, W2, NSigma;
    89   PhotCode *code;
    90 
    91   code = GetPhotcodebyCode (image[0].source);
    92114
    93115  /* reject multiple matched-stars */
     
    114136  ALLOCATE (dDmag, float, Ncal);
    115137  Nkeep = 0;
    116   Clam = code[0].C*0.001; /* photcode.C still in millimags */
    117   Klam = code[0].K;
    118   Xlam = code[0].X[0];
    119138  for (i = 0; i < Ncal; i++) {
    120139    /* if this entry has too many (or two few?) matches, skip it */
    121140    if (Nlist[Nstar[i]] != 1) continue;
    122     Mint = Clam + Mobs[i] + Xlam*Cref[i] + Klam*(Airm[i] - 1.000);
    123     Mabs = Mint - ZeroPt;
    124     /* note: subpix correction is applied in gstars */
    125141
    126     /* skip stars brighter than 8.0 */
    127     /* XXX EAM : perhaps make this more flexible?? */
    128     if (Mabs > 9.0) continue;
    129     if (Mabs < 5.0) continue;
     142    /* clip by instrumental magnitude */
     143    if (Minst[i] > CAL_INSTMAG_MAX) continue;
     144    if (Minst[i] < CAL_INSTMAG_MIN) continue;
    130145   
    131146    /* XXX EAM: note the artificial 0.005 dmag here */
     
    133148    dMo = MAX (0.005, dMobs[i]);
    134149
    135     Dmag[Nkeep] = (Mabs - Mref[i]);
     150    Dmag[Nkeep] = (Mobs[i] - Mref[i]);
    136151    dDmag[Nkeep] = (dMr*dMr + dMo*dMo);
    137152    Nkeep ++;
Note: See TracChangeset for help on using the changeset viewer.