Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/.mana =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/.mana (revision 37452) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/.mana (revision 37453) @@ -112,2 +112,28 @@ echo {$MAX / $dX^2} echo {3000 / 34} +echo {dsin(30+90)} {dcos(30)} +echo {dsin(30-90)} {-dcos(30)} +echo {dsin(30-90)} {-1*dcos(30)} +echo {dcos(30-90)} {dsin(30)} +echo {dcos(120-90)} {dsin(120)} +echo {dcos(-20-90)} {dsin(-20)} +echo {dcos(-20+90)} {dsin(-20)} +echo {dcos(30+90)} {dsin(30)} +echo {dcos(30+90)} {-1*dsin(30)} +echo {dcos(-20+90)} {-1*dsin(-20)} +!ps ux | grep fakes +!ps ux | grep fake +ll /proc/25713/fd/ +!ll /proc/25713/fd/ +!ll /proc/25713/fd/ +!ll /proc/25713/fd/ +!fhead /data/pikake.1/eugene/relastro.20140923/catdir.fake1/s3730.cpt +!ls -lrt /data/pikake.1/eugene/relastro.20140923/catdir.fake1/ +!ls -lrt /data/pikake.1/eugene/relastro.20140923/catdir.fake1/ +!ls -lrt /data/pikake.1/eugene/relastro.20140923/catdir.fake1/*.cpt +!ls -lrt /data/pikake.1/eugene/relastro.20140923/catdir.fake1/*.cpt +!ls -lrth /data/pikake.1/eugene/relastro.20140923/catdir.fake1/*.cpt +!ls -lrth /data/pikake.1/eugene/relastro.20140923/catdir.fake1/*.cp? +!ps ux | grep fake +kill 25713 +!kill 25713 Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/include/fakeastro.h =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/include/fakeastro.h (revision 37452) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/include/fakeastro.h (revision 37453) @@ -30,4 +30,5 @@ char CATDIR[DVO_MAX_PATH]; char SKY_TABLE[DVO_MAX_PATH]; +char MasterPhotcodeFile[DVO_MAX_PATH]; char CATMODE[16]; /* raw, mef, split, mysql */ @@ -46,4 +47,6 @@ int VERBOSE; int VERBOSE2; + +int FORCE; SkyRegion UserPatch; Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/ConfigInit.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/ConfigInit.c (revision 37452) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/ConfigInit.c (revision 37453) @@ -5,6 +5,4 @@ char *config, *file; char CatdirPhotcodeFile[256]; - struct stat filestat; - int status; /*** load configuration info ***/ @@ -30,4 +28,5 @@ ScanConfig(config, "CATMODE", "%s", 0, CATMODE); ScanConfig(config, "CATFORMAT", "%s", 0, CATFORMAT); + GetConfig (config, "PHOTCODE_FILE", "%s", 0, MasterPhotcodeFile); sprintf (ImageCat, "%s/Images.dat", CATDIR); @@ -39,14 +38,7 @@ if (*CATFORMAT == 0) strcpy (CATFORMAT, "ELIXIR"); - // check for existence of CATDIR - status = stat (CATDIR, &filestat); - if (status == -1) { - fprintf (stderr, "directory %s does not exist, giving up\n", CATDIR); - exit (1); - } - /* update master photcode table if not defined */ sprintf (CatdirPhotcodeFile, "%s/Photcodes.dat", CATDIR); - if (!LoadPhotcodes (CatdirPhotcodeFile, NULL, FALSE)) { + if (!LoadPhotcodes (CatdirPhotcodeFile, MasterPhotcodeFile, TRUE)) { fprintf (stderr, "error loading photcode table %s\n", CatdirPhotcodeFile); exit (1); Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/args.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/args.c (revision 37452) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/args.c (revision 37453) @@ -85,4 +85,10 @@ // } + FORCE = FALSE; + if ((N = get_argument (argc, argv, "-force"))) { + FORCE = TRUE; + remove_argument (N, &argc, argv); + } + VERBOSE = VERBOSE2 = FALSE; if ((N = get_argument (argc, argv, "-v"))) { Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_galaxy.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_galaxy.c (revision 37452) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_galaxy.c (revision 37453) @@ -4,6 +4,4 @@ int n, i; - - // INITTIME; // XXX check for an existing CATDIR and exit? @@ -18,13 +16,19 @@ gauss_init (50000); - int Nloop = 10; - int Nstars = 1000000; + int Nloop = 1; + int Nstars = 2000000; for (n = 0; n < Nloop; n++) { + INITTIME; + FakeAstro_Stars *stars = make_fakestars (Nstars); + + MARKTIME ("generate %d fake stars in %f sec\n", Nstars, dtime); // sort the stars by RA sortStars (stars, Nstars); + + MARKTIME ("sorted %d fake stars: %f sec\n", Nstars, dtime); // scan through the stars, loading the containing catalogs @@ -44,9 +48,17 @@ int Nsubset; FakeAstro_Stars *subset = make_subset (stars, Nstars, i, region, &Nsubset); + // MARKTIME ("subset %d fake stars: %f sec\n", Nsubset, dtime); // In parallel mode, write out the subset to a disk file. Block until a remote host // is available. In serial mode, just match against the appropriate region and save save_fakestars (subset, Nsubset, hosts, region, skylist[0].filename[0]); + // MARKTIME ("saved %d fake stars: %f sec\n", Nsubset, dtime); + + free (subset); + SkyListFree (skylist); } + free (stars); + + MARKTIME ("saved %d fake stars: %f sec\n", Nstars, dtime); // wait for last remote clients to finish Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_images.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_images.c (revision 37453) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_images.c (revision 37453) @@ -0,0 +1,37 @@ +# include "fakeastro.h" + +int fakeastro_images (SkyList *skylistInput, HostTable *hosts) { + + int Nimages; + Images *images = load_template_images (&Nimages); + + Images *bores = make_boresites (images, Nimages, Nbores); + + for (i = 0; i < Nbores; i++) { + + // we only want to make fake images for the exposures + if (strcmp(&images[i].coords.ctype[4], "-DIS")) continue; + + Images *fakeImages = make_fake_images (&images[i], &NfakeImages); + + for (j = 0; j < NfakeImages; j++) { + Stars *fakeStars = make_fake_stars (&fakeImages[i], &NfakeStars); + + // save fake stars... + } + } +} + +/* outline: + + * load reference image table + * generate a set of fake chips for each real exposure boresite + * select the real stars which land on the fake chips + * generate a set of measurements for the fake chip + * need to ingest the measurements into the new DVO as if they came from addstar + + * I can probably load stars for a decent fraction of the sky at once + */ + + +// -input catdir -output catdir -images images.fits Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/initialize.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/initialize.c (revision 37452) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/initialize.c (revision 37453) @@ -1,8 +1,26 @@ # include "fakeastro.h" + +void usage (void); +void usage_client (void); void initialize (int argc, char **argv) { + if (argc == 1) usage(); + if (get_argument (argc, argv, "-h")) usage(); + if (get_argument (argc, argv, "--h")) usage(); + if (get_argument (argc, argv, "-help")) usage(); + if (get_argument (argc, argv, "--help")) usage(); + ConfigInit (&argc, argv); args (argc, argv); + + // check for existence of CATDIR + struct stat filestat; + int status = stat (CATDIR, &filestat); + if (!FORCE && (status == 0)) { + fprintf (stderr, "directory %s exist, refusing to contaminate\n", CATDIR); + exit (1); + } + } Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/insert_fakestar.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/insert_fakestar.c (revision 37452) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/insert_fakestar.c (revision 37453) @@ -16,5 +16,7 @@ int Nsecfilt = catalog[0].Nsecfilt; - + NAVE = Nave = catalog[0].Naverage; + NSTARPAR = Nstarpar = catalog[0].Nstarpar; + /** incorporate unmatched image stars? **/ for (i = 0; i < Nstars; i++) { @@ -50,4 +52,7 @@ } + // I need a photcode for r-band + catalog[0].secfilt[Nave*Nsecfilt+1].M = stars[i].starpar.M_r + stars[i].starpar.DistMag; + catalog[0].starpar[Nstarpar] = stars[i].starpar; catalog[0].starpar[Nstarpar].averef = Nave; Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fake_images.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fake_images.c (revision 37453) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fake_images.c (revision 37453) @@ -0,0 +1,7 @@ +# include "fakeastro.h" + +Images *make_boresites (Images *images, int Nimages; int *nbores) { + + + +} Index: /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fakestars.c =================================================================== --- /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fakestars.c (revision 37452) +++ /branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fakestars.c (revision 37453) @@ -1,6 +1,14 @@ # include "fakeastro.h" -# define Z_GAL 50.0 +# define Z_GAL 500.0 # define R_GAL 2000.0 + +static double iFkap = 1.0 / 4.74047; // (km/sec/kpc) / (milli-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 FakeAstro_Stars *make_fakestars (int Nstars) { @@ -12,4 +20,5 @@ // but Greg Green / LSD use the older Reid et al 2004 (ApJ 616, 872) definition CoordTransform *transform = InitTransform (COORD_GALACTIC, COORD_CELESTIAL); + CoordTransform *transback = InitTransform (COORD_CELESTIAL, COORD_GALACTIC); int i; @@ -28,11 +37,20 @@ double z = rnd_gauss (0.0, Z_GAL); double r = sqrt(drand48()) * R_GAL; - double L = drand48() * 2 * M_PI; - - double distance = sqrt (SQ(r) + SQ(z)); - double B = atan2(z,r); + double Lrad = drand48() * 2 * M_PI; // double x = r*cos(L); // double y = r*sin(L); + + double distance = sqrt (SQ(r) + SQ(z)); + double Brad = atan2(z,r); + + double uL_gal = (A_oort * cos(2.0*Lrad) + B_oort) * cos(Brad) * iFkap; + double uB_gal = -0.5*A_oort * sin(2.0*Lrad) * sin(2.0*Brad) * iFkap; + + double uL_sol = (U_sol * sin(Lrad) - V_sol * cos(Lrad)) * iFkap / distance; + double uB_sol = ((U_sol * cos(Lrad) + V_sol * sin(Lrad))*sin(Brad) - W_sol*cos(Brad)) * iFkap / distance; + + double uL = uL_gal + uL_sol; + double uB = uB_gal + uB_sol; // crude Mr distribution from @@ -41,25 +59,48 @@ double Mr = rnd_gauss (11.25, 1.0); + double L = Lrad*DEG_RAD; + double B = Brad*DEG_RAD; + double R, D; ApplyTransform (&R, &D, L, B, transform); - stars[Nstars].R = R; - stars[Nstars].D = D; - stars[Nstars].flag = FALSE; - stars[Nstars].found = FALSE; + // C1, C2 are from http://arxiv.org/pdf/1306.2945v2.pdf + double Rrad = R*RAD_DEG; + double Drad = D*RAD_DEG; - stars[Nstars].starpar.galLon = L; - stars[Nstars].starpar.galLat = B; + double C1 = + cos(Drad)*transback->cos_phi + + sin(Drad)*cos(Rrad)*transback->sin_phi_sin_Xo - + sin(Drad)*sin(Rrad)*transback->sin_phi_cos_Xo; - stars[Nstars].starpar.Ebv = 0.0; - stars[Nstars].starpar.dEbv = 0.0; - stars[Nstars].starpar.DistMag = 5.0*log10(distance) - 5.0; - stars[Nstars].starpar.dDistMag = 0.0; - stars[Nstars].starpar.M_r = Mr; - stars[Nstars].starpar.dM_r = 0.0; - stars[Nstars].starpar.FeH = 0.0; - stars[Nstars].starpar.dFeH = 0.0; + double C2 = + - cos(Rrad)*transback->sin_phi_cos_Xo + - sin(Rrad)*transback->sin_phi_sin_Xo; - Nstars ++; + 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); + + stars[i].R = R; + stars[i].D = D; + stars[i].flag = FALSE; + stars[i].found = FALSE; + + stars[i].starpar.galLon = L; + stars[i].starpar.galLat = B; + + stars[i].starpar.Ebv = 0.0; + stars[i].starpar.dEbv = 0.0; + stars[i].starpar.DistMag = 5.0*log10(distance) - 5.0; + stars[i].starpar.dDistMag = 0.0; + stars[i].starpar.M_r = Mr; + stars[i].starpar.dM_r = 0.0; + stars[i].starpar.FeH = 0.0; + stars[i].starpar.dFeH = 0.0; + + stars[i].starpar.uRA = uR; + stars[i].starpar.uDEC = uD; } return stars;