IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 37602


Ignore:
Timestamp:
Nov 14, 2014, 8:40:19 PM (12 years ago)
Author:
eugene
Message:

move galaxy parameter functions to libdvo

Location:
branches/eam_branches/ipp-20140904/Ohana/src/uniphot
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20140904/Ohana/src/uniphot/include/setgalmodel.h

    r37434 r37602  
    1010char        *HOSTDIR;
    1111int          VERBOSE;
     12int          TESTING;
    1213int          UPDATE;
    1314int          PARALLEL;
    … …  
    3839// int           SetSignals                      PROTO((void));
    3940
    40 int           update_dvo_setgalmodel            PROTO((char *corrfile));
    41 int           update_dvo_setgalmodel_parallel   PROTO((SkyTable *sky, char *corrfile));
     41int           update_dvo_setgalmodel            PROTO((void));
     42int           update_dvo_setgalmodel_parallel   PROTO((SkyTable *sky));
    4243
    4344int           update_catalog_setgalmodel        PROTO((Catalog *catalog));

  • branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/initialize_setgalmodel.c

    r37434 r37602  
    22
    33void usage_setgalmodel () {
    4     fprintf (stderr, "USAGE: setgalmodel (corrections) [options]\n");
     4    fprintf (stderr, "USAGE: setgalmodel [options]\n");
    55    fprintf (stderr, "  options:\n");
    66    fprintf (stderr, "    -v : verbose mode\n");
    … …  
    4343    exit (1);
    4444  }
     45
     46  // XXX add to config?
     47  if (!InitGalaxyModel ("ROESER")) {
     48    fprintf (stderr, "failed to init galaxy model\n");
     49    exit (2);
     50  }
    4551}
    4652
    … …  
    8389  if ((N = get_argument (argc, argv, "-v"))) {
    8490    VERBOSE = TRUE;
     91    remove_argument (N, &argc, argv);
     92  }
     93
     94  TESTING = FALSE;
     95  if ((N = get_argument (argc, argv, "-testing"))) {
     96    TESTING = TRUE;
    8597    remove_argument (N, &argc, argv);
    8698  }
    … …  
    118130  }
    119131
    120   if (argc != 2) usage_setgalmodel();
     132  if (argc != 1) usage_setgalmodel();
    121133
    122134  return (TRUE);
    … …  
    124136
    125137void usage_setgalmodel_client () {
    126   fprintf (stderr, "USAGE: setgalmodel_client (correction) -hostID (hostID) -catdir (catdir) -hostdir (hostdir) [options]\n");
     138  fprintf (stderr, "USAGE: setgalmodel_client -hostID (hostID) -catdir (catdir) -hostdir (hostdir) [options]\n");
    127139  fprintf (stderr, "  options:\n");
    128140  fprintf (stderr, "    -region Rmin Rmax Dmin Dmax\n");
    … …  
    158170    exit (1);
    159171  }
     172
     173  // XXX add to config?
     174  if (!InitGalaxyModel ("ROESER")) {
     175    fprintf (stderr, "failed to init galaxy model\n");
     176    exit (2);
     177  }
    160178}
    161179
    … …  
    223241  if (!CATDIR) usage_setgalmodel_client();
    224242
    225   if (argc != 2) usage_setgalmodel_client();
     243  if (argc != 1) usage_setgalmodel_client();
    226244  return (TRUE);
    227245}

  • branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/setgalmodel.c

    r37434 r37602  
    88  initialize_setgalmodel (argc, argv);
    99
    10   // argv[1] holds the correction file
    11   status = update_dvo_setgalmodel (argv[1]);
     10  status = update_dvo_setgalmodel ();
    1211
    1312  if (!status) exit (1);

  • branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/setgalmodel_client.c

    r37434 r37602  
    1515  initialize_setgalmodel_client (argc, argv);
    1616
    17   update_dvo_setgalmodel (argv[1]);
     17  update_dvo_setgalmodel ();
    1818
    1919  exit (0);

  • branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/update_catalog_setgalmodel.c

    r37497 r37602  
    11# include "setgalmodel.h"
    2 
    3 static double iFkap = 0.001 / 4.74047; // (km/sec/kpc) / (arcsec/ year)
    4 static double A_oort = +14.5;  // km/sec/kpc
    5 static double B_oort = -13.0;  // km/sec/kpc
    6 
    7 static double U_sol  =  9.44; // km/sec
    8 static double V_sol  = 11.90; // km/sec
    9 static double W_sol  =  7.20; // km/sec
    102
    113int update_catalog_setgalmodel (Catalog *catalog) {
    … …  
    146
    157  // choose this number more judiciously...
    16   time_t timeRef = ohana_date_to_sec("2011/05/11,00:00:00");
     8  time_t timeRef = ohana_date_to_sec("2000/01/01,00:00:00");
    179
    1810  Average *average = catalog[0].average;
     11
     12  if ((catalog[0].Nstarpar_disk > 0) && (catalog[0].Nstarpar == 0)) {
     13    fprintf (stderr, "ERROR: failed to load starpar data\n");
     14    exit (2);
     15  }
    1916
    2017  CoordTransform *transform = InitTransform (COORD_CELESTIAL, COORD_GALACTIC);
    … …  
    2724    m = average[i].starparOffset;
    2825    StarPar *starpar = &catalog[0].starpar[m];
     26
     27    // fake or real QSOs are marked with FeH = +/- 100.0
     28    if (fabs(starpar->FeH) > 99.0) continue;
     29
    2930    double distance = pow(10.0, 0.2*(starpar->DistMag + 5.0));
    3031
    31     double glon = starpar->galLon * RAD_DEG;
    32     double glat = starpar->galLat * RAD_DEG;
     32    double Lrad = starpar->galLon * RAD_DEG;
     33    double Brad = starpar->galLat * RAD_DEG;
    3334
    34     double uL_gal = (A_oort * cos(2.0*glon) + B_oort) * cos(glat) * iFkap;
    35     double uB_gal = -0.5*A_oort * sin(2.0*glon) * sin(2.0*glat) * iFkap;
     35    double uL_gal, uB_gal;
     36    GalaxyMotionModel_radians(&uL_gal, &uB_gal, Lrad, Brad);
    3637
    37     double uL_sol =  (U_sol * sin(glon) - V_sol * cos(glon))                              * iFkap / distance;
    38     double uB_sol = ((U_sol * cos(glon) + V_sol * sin(glon))*sin(glat) - W_sol*cos(glat)) * iFkap / distance;
     38    double uL_sol, uB_sol;
     39    SolarMotionModel_radians(&uL_sol, &uB_sol, Lrad, Brad, distance);
    3940
    4041    double uL = uL_gal + uL_sol;
    4142    double uB = uB_gal + uB_sol;
    4243
    43     // C1, C2 are from http://arxiv.org/pdf/1306.2945v2.pdf
    44     double Rrad = average[i].R*RAD_DEG;
    45     double Drad = average[i].D*RAD_DEG;
    46 
    47     double C1 =
    48       cos(Drad)*transform->cos_phi +
    49       sin(Drad)*cos(Rrad)*transform->sin_phi_sin_Xo -
    50       sin(Drad)*sin(Rrad)*transform->sin_phi_cos_Xo;
    51 
    52     double C2 =
    53       - cos(Rrad)*transform->sin_phi_cos_Xo
    54       - sin(Rrad)*transform->sin_phi_sin_Xo;
    55 
    56     double cosBinv = 1.0 / sqrt(C1*C1 + C2*C2);
    57 
    58     // XXX add errors : I need to be able to choose the stars based on the error distribution
    59     double uR = cosBinv * (C1 * uL - C2 * uB);
    60     double uD = cosBinv * (C1 * uB + C2 * uL);
     44    // XXX: amplify motion to make tests easier:
     45    if (TESTING) {
     46      uL *= 100.0;
     47      uB *= 100.0;
     48    }
     49   
     50    double uR, uD;
     51    TransformProperMotion (&uR, &uD, uL, uB, average[i].R, average[i].D, transform);
    6152
    6253    starpar->uRA  = uR;
    … …  
    6960    for (j = 0; j < average[i].Nmeasure; j++) {
    7061
    71       int Nm = m + j;
    72       double dtime = (measure[Nm].t - timeRef) / 86400.0 / 365.25;
     62      double dtime = (measure[j].t - timeRef) / 86400.0 / 365.25;
    7363
    7464      // do not modify the original Xccd,Yccd values
    75       measure[Nm].RoffGAL = uR * dtime; // XXX divide by dcos? project and calculate?
    76       measure[Nm].DoffGAL = uD * dtime; // XXX divide by dcos? project and calculate?
     65      measure[j].RoffGAL = uR * dtime; // XXX divide by dcos? project and calculate?
     66      measure[j].DoffGAL = uD * dtime; // XXX divide by dcos? project and calculate?
    7767    }
    7868  }
    … …  
    8171}
    8272
    83 // I am using the galactic rotation and solar motion model to predict the reflex proper
    84 // motion of stars in the database.
    85 
    86 /*****
    87 
    88 From Siegfreid Roeser:
    89 
    90 Dear Gene,
    91 
    92  please find attached our approach for galactic rotation and
    93  solar motion that we assumed in our approach last fall.
    94 (Hope I understood correctly what you requested at the
    95 telecon). If not, please come back.
    96 
    97       fkap = 4.74047d-0   {km/s/kpc}/{mas/y}
    98 
    99       gl,gb  is galactic longitude resp. latitude
    100 
    101 
    102       PARAMETER (distance= 1.d0,Aoort=14.5d0 ,Boort=-13.d0  )
    103                       kpc           km/s/kpc   km/s/kpc
    104       PARAMETER (Usol = 9.44d0, Vsol =11.90d0, Wsol = 7.20d0)
    105                                    km/s
    106 C
    107 C      Gal. Rotation
    108 C
    109       proper motions in longitude/latitude
    110 
    111       emulgal = (AOORT*dcos(2.d0*gl)+BOORT)*dcos(gb)/fkap
    112       emubgal = -0.5d0*AOORT*dsin(2.d0*gl)*dsin(2.d0*gb)/fkap
    113 c
    114 c     Solar Motion
    115 cc
    116       emulsol = Usol*dsin(gl) - Vsol*dcos(gl)
    117       emulsol = emulsol/fkap/distance
    118       emubsol = Usol*dsin(gb)*dcos(gl)+Vsol*dsin(gb)*dsin(gl)
    119       emubsol = emubsol -Wsol*dcos(gb)
    120       emubsol = emubsol/fkap/distance
    121 
    122 
    123 ******/

  • branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/update_dvo_setgalmodel.c

    r37434 r37602  
    11# include "setgalmodel.h"
    22
    3 int update_dvo_setgalmodel (char *corrfile) {
     3int update_dvo_setgalmodel (void) {
    44
    55  SkyTable *sky = NULL;
    … …  
    1313 
    1414  if (PARALLEL && !HOST_ID) {
    15     int status = update_dvo_setgalmodel_parallel (sky, corrfile);
     15    int status = update_dvo_setgalmodel_parallel (sky);
    1616    return status;
    1717  }
    … …  
    3838    catalog.catformat = DVO_FORMAT_UNDEF; // not needed since we skip empty catalogs
    3939    catalog.catmode   = DVO_MODE_UNDEF;   // not needed since we skip empty catalogs
    40     catalog.catflags  = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
     40    catalog.catflags  = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF | LOAD_STARPAR;
    4141    catalog.Nsecfilt  = GetPhotcodeNsecfilt ();
    4242
    … …  
    7676# define DEBUG 1
    7777
    78 int update_dvo_setgalmodel_parallel (SkyTable *sky, char *corrfile) {
     78int update_dvo_setgalmodel_parallel (SkyTable *sky) {
    7979
    8080  // now launch the setgalmodel_client jobs to the parallel hosts
    81 
    82   char *abscorrfile = abspath (corrfile, DVO_MAX_PATH);
    8381
    8482  // load the list of hosts
    … …  
    9492
    9593    char command[1024];
    96     snprintf (command, 1024, "setgalmodel_client %s -hostID %d -catdir %s -hostdir %s -region %f %f %f %f",
    97               abscorrfile, table->hosts[i].hostID, CATDIR, table->hosts[i].pathname,
     94    snprintf (command, 1024, "setgalmodel_client -hostID %d -catdir %s -hostdir %s -region %f %f %f %f",
     95              table->hosts[i].hostID, CATDIR, table->hosts[i].pathname,
    9896              UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax);
    9997
    … …  
    124122    }
    125123  }
    126   free (abscorrfile);
    127124
    128125  // wait for the remote jobs to be completed
Note: See TracChangeset for help on using the changeset viewer.