Index: /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/include/setgalmodel.h =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/include/setgalmodel.h (revision 37601) +++ /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/include/setgalmodel.h (revision 37602) @@ -10,4 +10,5 @@ char *HOSTDIR; int VERBOSE; +int TESTING; int UPDATE; int PARALLEL; @@ -38,6 +39,6 @@ // int SetSignals PROTO((void)); -int update_dvo_setgalmodel PROTO((char *corrfile)); -int update_dvo_setgalmodel_parallel PROTO((SkyTable *sky, char *corrfile)); +int update_dvo_setgalmodel PROTO((void)); +int update_dvo_setgalmodel_parallel PROTO((SkyTable *sky)); int update_catalog_setgalmodel PROTO((Catalog *catalog)); Index: /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/initialize_setgalmodel.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/initialize_setgalmodel.c (revision 37601) +++ /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/initialize_setgalmodel.c (revision 37602) @@ -2,5 +2,5 @@ void usage_setgalmodel () { - fprintf (stderr, "USAGE: setgalmodel (corrections) [options]\n"); + fprintf (stderr, "USAGE: setgalmodel [options]\n"); fprintf (stderr, " options:\n"); fprintf (stderr, " -v : verbose mode\n"); @@ -43,4 +43,10 @@ exit (1); } + + // XXX add to config? + if (!InitGalaxyModel ("ROESER")) { + fprintf (stderr, "failed to init galaxy model\n"); + exit (2); + } } @@ -83,4 +89,10 @@ if ((N = get_argument (argc, argv, "-v"))) { VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + + TESTING = FALSE; + if ((N = get_argument (argc, argv, "-testing"))) { + TESTING = TRUE; remove_argument (N, &argc, argv); } @@ -118,5 +130,5 @@ } - if (argc != 2) usage_setgalmodel(); + if (argc != 1) usage_setgalmodel(); return (TRUE); @@ -124,5 +136,5 @@ void usage_setgalmodel_client () { - fprintf (stderr, "USAGE: setgalmodel_client (correction) -hostID (hostID) -catdir (catdir) -hostdir (hostdir) [options]\n"); + fprintf (stderr, "USAGE: setgalmodel_client -hostID (hostID) -catdir (catdir) -hostdir (hostdir) [options]\n"); fprintf (stderr, " options:\n"); fprintf (stderr, " -region Rmin Rmax Dmin Dmax\n"); @@ -158,4 +170,10 @@ exit (1); } + + // XXX add to config? + if (!InitGalaxyModel ("ROESER")) { + fprintf (stderr, "failed to init galaxy model\n"); + exit (2); + } } @@ -223,5 +241,5 @@ if (!CATDIR) usage_setgalmodel_client(); - if (argc != 2) usage_setgalmodel_client(); + if (argc != 1) usage_setgalmodel_client(); return (TRUE); } Index: /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/setgalmodel.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/setgalmodel.c (revision 37601) +++ /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/setgalmodel.c (revision 37602) @@ -8,6 +8,5 @@ initialize_setgalmodel (argc, argv); - // argv[1] holds the correction file - status = update_dvo_setgalmodel (argv[1]); + status = update_dvo_setgalmodel (); if (!status) exit (1); Index: /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/setgalmodel_client.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/setgalmodel_client.c (revision 37601) +++ /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/setgalmodel_client.c (revision 37602) @@ -15,5 +15,5 @@ initialize_setgalmodel_client (argc, argv); - update_dvo_setgalmodel (argv[1]); + update_dvo_setgalmodel (); exit (0); Index: /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/update_catalog_setgalmodel.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/update_catalog_setgalmodel.c (revision 37601) +++ /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/update_catalog_setgalmodel.c (revision 37602) @@ -1,11 +1,3 @@ # include "setgalmodel.h" - -static double iFkap = 0.001 / 4.74047; // (km/sec/kpc) / (arcsec/ year) -static double A_oort = +14.5; // km/sec/kpc -static double B_oort = -13.0; // km/sec/kpc - -static double U_sol = 9.44; // km/sec -static double V_sol = 11.90; // km/sec -static double W_sol = 7.20; // km/sec int update_catalog_setgalmodel (Catalog *catalog) { @@ -14,7 +6,12 @@ // choose this number more judiciously... - time_t timeRef = ohana_date_to_sec("2011/05/11,00:00:00"); + time_t timeRef = ohana_date_to_sec("2000/01/01,00:00:00"); Average *average = catalog[0].average; + + if ((catalog[0].Nstarpar_disk > 0) && (catalog[0].Nstarpar == 0)) { + fprintf (stderr, "ERROR: failed to load starpar data\n"); + exit (2); + } CoordTransform *transform = InitTransform (COORD_CELESTIAL, COORD_GALACTIC); @@ -27,36 +24,30 @@ m = average[i].starparOffset; StarPar *starpar = &catalog[0].starpar[m]; + + // fake or real QSOs are marked with FeH = +/- 100.0 + if (fabs(starpar->FeH) > 99.0) continue; + double distance = pow(10.0, 0.2*(starpar->DistMag + 5.0)); - double glon = starpar->galLon * RAD_DEG; - double glat = starpar->galLat * RAD_DEG; + double Lrad = starpar->galLon * RAD_DEG; + double Brad = starpar->galLat * RAD_DEG; - double uL_gal = (A_oort * cos(2.0*glon) + B_oort) * cos(glat) * iFkap; - double uB_gal = -0.5*A_oort * sin(2.0*glon) * sin(2.0*glat) * iFkap; + double uL_gal, uB_gal; + GalaxyMotionModel_radians(&uL_gal, &uB_gal, Lrad, Brad); - double uL_sol = (U_sol * sin(glon) - V_sol * cos(glon)) * iFkap / distance; - double uB_sol = ((U_sol * cos(glon) + V_sol * sin(glon))*sin(glat) - W_sol*cos(glat)) * iFkap / distance; + double uL_sol, uB_sol; + SolarMotionModel_radians(&uL_sol, &uB_sol, Lrad, Brad, distance); double uL = uL_gal + uL_sol; double uB = uB_gal + uB_sol; - // C1, C2 are from http://arxiv.org/pdf/1306.2945v2.pdf - double Rrad = average[i].R*RAD_DEG; - double Drad = average[i].D*RAD_DEG; - - double C1 = - cos(Drad)*transform->cos_phi + - sin(Drad)*cos(Rrad)*transform->sin_phi_sin_Xo - - sin(Drad)*sin(Rrad)*transform->sin_phi_cos_Xo; - - double C2 = - - cos(Rrad)*transform->sin_phi_cos_Xo - - sin(Rrad)*transform->sin_phi_sin_Xo; - - double cosBinv = 1.0 / sqrt(C1*C1 + C2*C2); - - // XXX add errors : I need to be able to choose the stars based on the error distribution - double uR = cosBinv * (C1 * uL - C2 * uB); - double uD = cosBinv * (C1 * uB + C2 * uL); + // XXX: amplify motion to make tests easier: + if (TESTING) { + uL *= 100.0; + uB *= 100.0; + } + + double uR, uD; + TransformProperMotion (&uR, &uD, uL, uB, average[i].R, average[i].D, transform); starpar->uRA = uR; @@ -69,10 +60,9 @@ for (j = 0; j < average[i].Nmeasure; j++) { - int Nm = m + j; - double dtime = (measure[Nm].t - timeRef) / 86400.0 / 365.25; + double dtime = (measure[j].t - timeRef) / 86400.0 / 365.25; // do not modify the original Xccd,Yccd values - measure[Nm].RoffGAL = uR * dtime; // XXX divide by dcos? project and calculate? - measure[Nm].DoffGAL = uD * dtime; // XXX divide by dcos? project and calculate? + measure[j].RoffGAL = uR * dtime; // XXX divide by dcos? project and calculate? + measure[j].DoffGAL = uD * dtime; // XXX divide by dcos? project and calculate? } } @@ -81,43 +71,2 @@ } -// I am using the galactic rotation and solar motion model to predict the reflex proper -// motion of stars in the database. - -/***** - -From Siegfreid Roeser: - -Dear Gene, - - please find attached our approach for galactic rotation and - solar motion that we assumed in our approach last fall. -(Hope I understood correctly what you requested at the -telecon). If not, please come back. - - fkap = 4.74047d-0 {km/s/kpc}/{mas/y} - - gl,gb is galactic longitude resp. latitude - - - PARAMETER (distance= 1.d0,Aoort=14.5d0 ,Boort=-13.d0 ) - kpc km/s/kpc km/s/kpc - PARAMETER (Usol = 9.44d0, Vsol =11.90d0, Wsol = 7.20d0) - km/s -C -C Gal. Rotation -C - proper motions in longitude/latitude - - emulgal = (AOORT*dcos(2.d0*gl)+BOORT)*dcos(gb)/fkap - emubgal = -0.5d0*AOORT*dsin(2.d0*gl)*dsin(2.d0*gb)/fkap -c -c Solar Motion -cc - emulsol = Usol*dsin(gl) - Vsol*dcos(gl) - emulsol = emulsol/fkap/distance - emubsol = Usol*dsin(gb)*dcos(gl)+Vsol*dsin(gb)*dsin(gl) - emubsol = emubsol -Wsol*dcos(gb) - emubsol = emubsol/fkap/distance - - -******/ Index: /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/update_dvo_setgalmodel.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/update_dvo_setgalmodel.c (revision 37601) +++ /branches/eam_branches/ipp-20140904/Ohana/src/uniphot/src/update_dvo_setgalmodel.c (revision 37602) @@ -1,5 +1,5 @@ # include "setgalmodel.h" -int update_dvo_setgalmodel (char *corrfile) { +int update_dvo_setgalmodel (void) { SkyTable *sky = NULL; @@ -13,5 +13,5 @@ if (PARALLEL && !HOST_ID) { - int status = update_dvo_setgalmodel_parallel (sky, corrfile); + int status = update_dvo_setgalmodel_parallel (sky); return status; } @@ -38,5 +38,5 @@ catalog.catformat = DVO_FORMAT_UNDEF; // not needed since we skip empty catalogs catalog.catmode = DVO_MODE_UNDEF; // not needed since we skip empty catalogs - catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF; + catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF | LOAD_STARPAR; catalog.Nsecfilt = GetPhotcodeNsecfilt (); @@ -76,9 +76,7 @@ # define DEBUG 1 -int update_dvo_setgalmodel_parallel (SkyTable *sky, char *corrfile) { +int update_dvo_setgalmodel_parallel (SkyTable *sky) { // now launch the setgalmodel_client jobs to the parallel hosts - - char *abscorrfile = abspath (corrfile, DVO_MAX_PATH); // load the list of hosts @@ -94,6 +92,6 @@ char command[1024]; - snprintf (command, 1024, "setgalmodel_client %s -hostID %d -catdir %s -hostdir %s -region %f %f %f %f", - abscorrfile, table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, + snprintf (command, 1024, "setgalmodel_client -hostID %d -catdir %s -hostdir %s -region %f %f %f %f", + table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax); @@ -124,5 +122,4 @@ } } - free (abscorrfile); // wait for the remote jobs to be completed