IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 37564


Ignore:
Timestamp:
Nov 6, 2014, 8:18:39 PM (12 years ago)
Author:
eugene
Message:

added 2mass, more testing

Location:
branches/eam_branches/ipp-20140904/Ohana/src/fakeastro
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/include/fakeastro.h

    r37532 r37564  
    8787int    VERBOSE;
    8888int    VERBOSE2;
     89int    TESTING;
    8990
    9091int    FORCE;
     92int    UNIFORM_RADEC;
    9193
    92 int FAKEASTRO_NLOOP;
    93 int FAKEASTRO_NSTARS;
    94 int FAKEASTRO_NQSO_ICRF;
    95 int FAKEASTRO_NQSO_ZERO;
    96 float FAKEASTRO_ZGAL; // parsec
    97 float FAKEASTRO_RGAL; // parsec
    98 char FAKEASTRO_REF_EPOCH[80];
    99 char FAKEASTRO_2MASS_EPOCH[80];
     94int    FAKEASTRO_NLOOP;
     95int    FAKEASTRO_NSTARS;
     96int    FAKEASTRO_NQSO_ICRF;
     97int    FAKEASTRO_NQSO_ZERO;
     98float  FAKEASTRO_ZGAL; // parsec
     99float  FAKEASTRO_RGAL; // parsec
     100char   FAKEASTRO_REF_EPOCH[80];
     101char   FAKEASTRO_2MASS_EPOCH[80];
    100102
    101103float  RADIUS;
  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/ConfigInit.c

    r37523 r37564  
    1919  // if (!ScanConfig (config, "ADDSTAR_RADIUS",         "%lf", 0, &ADDSTAR_RADIUS))   ADDSTAR_RADIUS = 1.0;
    2020
    21   if (FAKEASTRO_OP == OP_GALAXY) {
     21  if ((FAKEASTRO_OP == OP_GALAXY) || (FAKEASTRO_OP == OP_2MASS)) {
    2222    // force CATDIR to be absolute (so parallel mode will work)
    2323    GetConfig (config, "CATDIR",                 "%s",  0, CATDIR);
     
    7373  if (*CATFORMAT == 0) strcpy (CATFORMAT, "ELIXIR");
    7474
    75   char *CATDIR_CHECK = (FAKEASTRO_OP == OP_GALAXY) ? CATDIR : CATDIR_OUTPUT;
     75  char *CATDIR_CHECK = (FAKEASTRO_OP == OP_IMAGES) ? CATDIR_OUTPUT : CATDIR;
    7676
    77   // check for existence of CATDIR
    78   struct stat filestat;
    79   int status = stat (CATDIR_CHECK, &filestat);
    80   if (!FORCE && (status == 0)) {
    81     fprintf (stderr, "directory %s exists, refusing to contaminate\n", CATDIR);
    82     exit (1);
     77  // OP_2MASS is adding detections to an existing db, the others require and empty db
     78  if (FAKEASTRO_OP != OP_2MASS) {
     79    // check for existence of CATDIR
     80    struct stat filestat;
     81    int status = stat (CATDIR_CHECK, &filestat);
     82    if (!FORCE && (status == 0)) {
     83      fprintf (stderr, "directory %s exists, refusing to contaminate\n", CATDIR);
     84      exit (1);
     85    }
    8386  }
    8487
  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/args.c

    r37463 r37564  
    1919    remove_argument (N, argc, argv);
    2020    FAKEASTRO_OP = OP_GALAXY;
     21  }
     22
     23  if ((N = get_argument (*argc, argv, "-2mass"))) {
     24    remove_argument (N, argc, argv);
     25    FAKEASTRO_OP = OP_2MASS;
     26  }
     27
     28  TESTING = FALSE;
     29  if ((N = get_argument (*argc, argv, "-testing"))) {
     30    remove_argument (N, argc, argv);
     31    TESTING = TRUE;
    2132  }
    2233
     
    8495  //   MaxDensityUse = TRUE;
    8596  // }
     97
     98  UNIFORM_RADEC = FALSE;
     99  if ((N = get_argument (*argc, argv, "-uniform-radec"))) {
     100    UNIFORM_RADEC = TRUE;
     101    remove_argument (N, argc, argv);
     102  }
    86103
    87104  FORCE = FALSE;
  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_2mass.c

    r37523 r37564  
    2323    catalog.catflags  = LOAD_AVES | LOAD_SECF | LOAD_MEAS | LOAD_STARPAR;
    2424    catalog.Nsecfilt  = GetPhotcodeNsecfilt ();
    25     if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "r")) {
     25    if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "w")) {
    2626      fprintf (stderr, "ERROR: failure reading catalog %s\n", catalog.filename);
    2727      exit (1);
     
    3737
    3838    SetProtect (TRUE);
    39     dvo_catalog_update (&catalog, VERBOSE);
     39    if (!dvo_catalog_update (&catalog, VERBOSE)) {
     40      fprintf (stderr, "ERROR: failure to update %s\n", catalog.filename);
     41      exit (3);
     42    }
    4043    SetProtect (FALSE);
    4144
    42     dvo_catalog_unlock (&catalog);
     45    if (!dvo_catalog_unlock (&catalog)) {
     46      fprintf (stderr, "ERROR: failure to unlock %s\n", catalog.filename);
     47      exit (2);
     48    }
    4349    dvo_catalog_free (&catalog);
    4450    if (VERBOSE) MARKTIME ("save cpt: %f sec\n", dtime); RESETTIME;
  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_2mass_measures.c

    r37523 r37564  
    4444    double J_PS1 = secfilt[i*Nsecfilt + Nsec].M - 1.0; // make all stars have y-J = 1.0
    4545    if (J_PS1 > 16.0) continue; // only generate 2MASS detections for objects with J_PS1 < 16.0
     46
     47    if (isnan(secfilt[i*Nsecfilt + Nsec].M)) {
     48      // look for a non-NAN secfilt mag and just use that (it is not super important)
     49      int ns;
     50      for (ns = 0; ns < Nsecfilt; ns++) {
     51        if (!isnan(secfilt[i*Nsecfilt + ns].M)) break;
     52      }
     53      if (ns == Nsecfilt) continue; // no non-nan
     54
     55      J_PS1 = secfilt[i*Nsecfilt + ns].M; // pretend secfilt.M = J
     56      if (J_PS1 > 16.0) continue; // only generate 2MASS detections for objects with J_PS1 < 16.0
     57    }
    4658
    4759    double Minst = J_PS1 - ZP;
     
    101113
    102114    // This is may optionally be replaced by the internal sequence (see FilterStars.c)
    103     measure[Nmeasure].detID      = Nmeasure;
     115    measure[Nmeasure].detID      = 0;
    104116    Nmeasure ++;
    105117  }
  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fakeqsos.c

    r37527 r37564  
    1212  float dPmax = sin(RAD_DEG*UserPatch.Dmax);
    1313  float dP = dPmax - dPmin;
     14
     15  CoordTransform *transback = InitTransform (COORD_CELESTIAL, COORD_GALACTIC);
    1416
    1517  int i;
     
    3234      break;
    3335    }
     36    double L, B;
     37    ApplyTransform (&L, &B, R, D, transback);
    3438
    3539    // Mr will be interpretted as m_r
    36     double Mr = rnd_gauss (17.0, 1.0);
     40    double Mr = rnd_gauss (18.0, 1.5);
    3741   
    3842    stars[i].R = R;
     
    4347    stars[i].starpar.R      = R;
    4448    stars[i].starpar.D      = D;
    45     stars[i].starpar.galLon = NAN;
    46     stars[i].starpar.galLat = NAN;
     49    stars[i].starpar.galLon = L;
     50    stars[i].starpar.galLat = B;
    4751
    4852    // how shall I distinguish ICRF and ZERO quasars
  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fakestars.c

    r37523 r37564  
    1919  CoordTransform *transback = InitTransform (COORD_CELESTIAL, COORD_GALACTIC);
    2020
     21  float dR = UserPatch.Rmax - UserPatch.Rmin;
     22  float dPmin = sin(RAD_DEG*UserPatch.Dmin);
     23  float dPmax = sin(RAD_DEG*UserPatch.Dmax);
     24  float dP = dPmax - dPmin;
     25
    2126  int i;
    2227  for (i = 0; i < Nstars; i++) {
     
    3540    double z,r,L,B,R,D,Lrad,Brad;
    3641
    37     int inPatch = FALSE;
    38     while (!inPatch) {
    39       z = rnd_gauss (0.0, FAKEASTRO_ZGAL);
    40       r = sqrt(drand48()) * FAKEASTRO_RGAL;
    41       Lrad = drand48() * 2 * M_PI;
    42       Brad = atan2(z,r);
     42    if (UNIFORM_RADEC) {
     43      // we can generate a distribution which is uniform on the sky, in which case
     44      // we can limit to the selected patch analyically
     45      r = pow(drand48(), 0.33333) * FAKEASTRO_RGAL;
     46      z = 0.0;
     47      double phi = dP * drand48() + dPmin;
     48      D = DEG_RAD * asin(phi); // random in degrees
     49      R = drand48() * dR + UserPatch.Rmin;   // random in degrees
     50      R = ohana_normalize_angle (R);
     51
     52      ApplyTransform (&L, &B, R, D, transback);
     53      Lrad = L*RAD_DEG;
     54      Brad = B*RAD_DEG;
     55    } else {
     56      int inPatch = FALSE;
     57      while (!inPatch) {
     58        z = rnd_gauss (0.0, FAKEASTRO_ZGAL);
     59        r = sqrt(drand48()) * FAKEASTRO_RGAL;
     60        Lrad = drand48() * 2 * M_PI;
     61        Brad = atan2(z,r);
    4362     
    44       L = Lrad*DEG_RAD;
    45       B = Brad*DEG_RAD;
     63        L = Lrad*DEG_RAD;
     64        B = Brad*DEG_RAD;
    4665     
    47       ApplyTransform (&R, &D, L, B, transform);
    48       if (R < UserPatch.Rmin) continue;
    49       if (R > UserPatch.Rmax) continue;
    50       if (D < UserPatch.Dmin) continue;
    51       if (D > UserPatch.Dmax) continue;
    52       break;
     66        ApplyTransform (&R, &D, L, B, transform);
     67        if (R < UserPatch.Rmin) continue;
     68        if (R > UserPatch.Rmax) continue;
     69        if (D < UserPatch.Dmin) continue;
     70        if (D > UserPatch.Dmax) continue;
     71        break;
     72      }
    5373    }
    5474
     
    6686    double uL = uL_gal + uL_sol;
    6787    double uB = uB_gal + uB_sol;
    68 
     88   
     89    // XXX: amplify motion to make tests easier:
     90    if (TESTING) {
     91      uL *= 100.0;
     92      uB *= 100.0;
     93    }
     94   
    6995    // crude Mr distribution from Bochanski et al 2010
    7096    // http://iopscience.iop.org/1538-3881/139/6/2679/pdf/aj_139_6_2679.pdf
    71 
    72     double Mr = rnd_gauss (11.25, 1.0);
     97   
     98    // two gaussian distributions:
     99    // 75% in narrow gauss, with Mr = 11.25, sigma = 1.0
     100    // 25% is wide gauss, with Mr = 10.0, sigma = 3.0
     101    // first choose which gauss:
     102    int bigPeak = (drand48() > 0.60);
     103    double Mr;
     104    if (bigPeak) {
     105      Mr = rnd_gauss (11.25, 1.0);
     106    } else {
     107      Mr = rnd_gauss (10.00, 2.5);
     108    }
    73109   
    74110    // C1, C2 are from http://arxiv.org/pdf/1306.2945v2.pdf
  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/match_fake_stars.c

    r37532 r37564  
    252252    dvo_starpar_init (&catalog[0].starpar[Nstarpar]);
    253253    catalog[0].starpar[Nstarpar] = stars[i].starpar;
     254    catalog[0].starpar[Nstarpar].averef   = Nave;
     255    catalog[0].starpar[Nstarpar].objID    = catalog[0].average[Nave].objID;
     256    catalog[0].starpar[Nstarpar].catID    = catalog[0].catID;
    254257    Nstarpar ++;
    255258
Note: See TracChangeset for help on using the changeset viewer.