Changeset 37422

Timestamp:

Sep 23, 2014, 9:53:58 AM (12 years ago)

Author:

eugene

Message:

loadstarpar basically complete, not tested

Location:

branches/eam_branches/ipp-20140904/Ohana/src/addstar

Files:

• Makefile (modified) (2 diffs)
• include/loadstarpar.h (modified) (1 diff)
• src/args_loadstarpar.c (modified) (1 diff)
• src/coord_transforms.c (modified) (2 diffs)
• src/find_matches_starpar.c (modified) (2 diffs)
• src/loadstarpar.c (modified) (1 diff)
• src/loadstarpar_make_subset.c (modified) (1 diff)
• src/loadstarpar_readstars.c (modified) (5 diffs)
• src/loadstarpar_save_remote.c (modified) (2 diffs)

branches/eam_branches/ipp-20140904/Ohana/src/addstar/Makefile

@@ -250 +250 @@
 $(SRC)/SkyRegionUtils.$(ARCH).o \
 $(SRC)/args_loadstarpar.$(ARCH).o \
-$(SRC)/coord_transforms.$(ARCH).o \
 $(SRC)/find_matches_starpar.$(ARCH).o \
 $(SRC)/loadstarpar_catalog.$(ARCH).o \
@@ -268 +267 @@
 $(SRC)/SkyRegionUtils.$(ARCH).o \
 $(SRC)/args_loadstarpar.$(ARCH).o \
-$(SRC)/coord_transforms.$(ARCH).o \
 $(SRC)/find_matches_starpar.$(ARCH).o \
 $(SRC)/loadstarpar_catalog.$(ARCH).o \

branches/eam_branches/ipp-20140904/Ohana/src/addstar/include/loadstarpar.h

@@ -3 +3 @@
   double R, D;
   StarPar starpar;
-  int flag; // XXX how is this used?
+  int flag; // in a subset?
   int found; // assigned to an object?
 } StarPar_Stars;

branches/eam_branches/ipp-20140904/Ohana/src/addstar/src/args_loadstarpar.c

@@ -23 +23 @@
   UserPatch.Dmin = -90;
   UserPatch.Dmax = +90;
+
+  // XXX for the moment, make this selection manual.  it needs to be automatic 
+  // based on the state of the SkyTable
+  HOST_ID = 0;
+  PARALLEL = FALSE;
+  if ((N = get_argument (argc, argv, "-parallel"))) {
+    PARALLEL = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+  // this is a test mode : rather than launching the remote jobs and waiting for completion,
+  // relphot will simply list the remote command and wait for the user to signal completion
+  PARALLEL_MANUAL = FALSE;
+  if ((N = get_argument (argc, argv, "-parallel-manual"))) {
+    PARALLEL = TRUE; // -parallel-manual implies -parallel
+    PARALLEL_MANUAL = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+  // this is a test mode : rather than launching the relphot_client jobs remotely, they are 
+  // run in serial via 'system'
+  PARALLEL_SERIAL = FALSE;
+  if ((N = get_argument (argc, argv, "-parallel-serial"))) {
+    if (PARALLEL_MANUAL) {
+      fprintf (stderr, "ERROR: cannot mix -parallel-manual and -parallel-serial\n");
+      exit (1);
+    }
+    PARALLEL = TRUE; // -parallel-serial implies -parallel
+    PARALLEL_SERIAL = TRUE;
+    remove_argument (N, &argc, argv);
+  }
 
   /* only add to existing regions */

branches/eam_branches/ipp-20140904/Ohana/src/addstar/src/coord_transforms.c

@@ -1 +1 @@
 # include "addstar.h"
 # include "loadstarpar.h"
+
+// values used by libdvo/coord_systems (from Liu et al 2011, A&A 526, A16):
+// Liu etal and Reid etal disagree at the 0.3 arcsec level
+// transform->phi =  62.8717488056*RAD_DEG;
+// transform->Xo  =  32.9319185700*RAD_DEG;
+// transform->xo  = 282.8594812080;
 
 int galactic_to_celestial (double *R, double *D, double l, double b) {
@@ -7 +13 @@
   abort();
 }
+
+int celestial_to_galactic (double *l, double *b, double R, double D) {
+
+_angp = np.radians(192.859508333) #  12h 51m 26.282s (J2000)
+_dngp = np.radians(27.128336111)  # +27d 07' 42.01" (J2000) 
+_l0   = np.radians(32.932)
+_ce   = np.cos(_dngp)
+_se   = np.sin(_dngp)
+
+
+
+}
+
+
+
+
+/*** python version from LSD (Mario Juric)
+
+# Appendix of Reid et al. (http://adsabs.harvard.edu/cgi-bin/bib_query?2004ApJ...616..872R)
+# This convention is also used by LAMBDA/WMAP (http://lambda.gsfc.nasa.gov/toolbox/tb_coordconv.cfm)
+_angp = np.radians(192.859508333) #  12h 51m 26.282s (J2000)
+_dngp = np.radians(27.128336111)  # +27d 07' 42.01" (J2000) 
+_l0   = np.radians(32.932)
+_ce   = np.cos(_dngp)
+_se   = np.sin(_dngp)
+
+def equgal(ra, dec):
+        ra = np.radians(ra)
+        dec = np.radians(dec)
+
+        cd, sd = np.cos(dec), np.sin(dec)
+        ca, sa = np.cos(ra - _angp), np.sin(ra - _angp)
+
+        sb = cd*_ce*ca + sd*_se
+        l = np.arctan2(sd - sb*_se, cd*sa*_ce) + _l0
+        b = np.arcsin(sb)
+
+        l = np.where(l < 0, l + 2.*np.pi, l)
+
+        l = np.degrees(l)
+        b = np.degrees(b)
+
+        return NamedList(('l', l), ('b', b))
+
+def galequ(l, b):
+        l = np.radians(l)
+        b = np.radians(b)
+
+        cb, sb = np.cos(b), np.sin(b)
+        cl, sl = np.cos(l-_l0), np.sin(l-_l0)
+
+        ra  = np.arctan2(cb*cl, sb*_ce-cb*_se*sl) + _angp
+        dec = np.arcsin(cb*_ce*sl + sb*_se)
+
+        ra = np.where(ra < 0, ra + 2.*np.pi, ra)
+
+        ra = np.degrees(ra)
+        dec = np.degrees(dec)
+
+        return NamedList(('ra', ra), ('dec', dec))
+
+ ***/

branches/eam_branches/ipp-20140904/Ohana/src/addstar/src/find_matches_starpar.c

@@ -151 +151 @@
   }
 
-# if (0)
+# if (1)
   /** incorporate unmatched image stars? **/
   for (i = 0; (i < Nstars) && !options->only_match; i++) {
@@ -157 +157 @@
     // skip already matched stars
     if (stars[i].found != -1) continue;
-    if (!IN_REGION (stars[i].average.R, stars[i].average.D)) continue;
+    if (!IN_REGION (stars[i].R, stars[i].D)) continue;
 
     /* make sure there is space for next entry */
     if (Nstarpar >= NSTARPAR) {
       NSTARPAR = Nstarpar + 1000;
-      REALLOCATE (catalog[0].starpar, Starpar, NSTARPAR);
+      REALLOCATE (catalog[0].starpar, StarPar, NSTARPAR);
     }
     if (Nave >= NAVE) {

branches/eam_branches/ipp-20140904/Ohana/src/addstar/src/loadstarpar.c

@@ -25 +25 @@
   
   // load the list of hosts
-  HostTable *hosts = HostTableLoad (CATDIR, sky->hosts);
-  if (!hosts) {
-    fprintf (stderr, "ERROR: failure reading Host Table %s for database %s\n", sky->hosts, CATDIR);
-    exit (1);
-  }    
+  HostTable *hosts = NULL;
+  if (PARALLEL) {
+    hosts = HostTableLoad (CATDIR, sky->hosts);
+    if (!hosts) {
+      fprintf (stderr, "ERROR: failure reading Host Table %s for database %s\n", sky->hosts, CATDIR);
+      exit (1);
+    }    
+
+    // ensure that the paths are absolute path names
+    for (i = 0; i < hosts->Nhosts; i++) {
+      char *tmppath = abspath (hosts->hosts[i].pathname, DVO_MAX_PATH);
+      free (hosts->hosts[i].pathname);
+      hosts->hosts[i].pathname = tmppath;
+    }
+
+    // set up the array of active hosts
+    init_remote_hosts ();
+  }
 
   // generate the subset matching the user-selected region

branches/eam_branches/ipp-20140904/Ohana/src/addstar/src/loadstarpar_make_subset.c

@@ -35 +35 @@
     Nsubset ++;
 
+    stars[i].flag = TRUE;
+
     CHECK_REALLOCATE (subset, StarPar_Stars, NSUBSET, Nsubset, 10000);
   }

branches/eam_branches/ipp-20140904/Ohana/src/addstar/src/loadstarpar_readstars.c

@@ -44 +44 @@
     return (NULL);
   }
-  fclose (f);
 
   char type[16];
@@ -75 +74 @@
     return NULL;
   }
+  fclose (f);
 
   double *errImage = (double *) matrix.buffer;
@@ -93 +93 @@
   ALLOCATE (stars, StarPar_Stars, NSTARS);
 
+  // libdvo uses Liu et al 2011 (A&A 526, A16) for default galactic coords,
+  // but Greg Green / LSD use the older Reid et al 2004 (ApJ 616, 872) definition
+  CoordTransform *transform = InitTransform (COORD_GALACTIC_REID_2004, COORD_CELESTIAL);
+
   for (i = 0; i < NstarsIn; i++) {
 
@@ -100 +104 @@
 
     double R, D;
-    galactic_to_celestial (&R, &D, glon[i], glat[i]);
+    ApplyTransform (&R, &D, glon[i], glat[i], transform);
 
     Rmin = MIN (Rmin, R);
@@ -114 +118 @@
     stars[Nstars].R = R;
     stars[Nstars].D = D;
+    stars[Nstars].flag  = FALSE;
+    stars[Nstars].found = FALSE;
 
     stars[Nstars].starpar.galLon = glon[i];

branches/eam_branches/ipp-20140904/Ohana/src/addstar/src/loadstarpar_save_remote.c

@@ -24 +24 @@
 
   // if this region is a parallel thing, save and launch remote
-  if (!region->hostID) { 
+  if (!PARALLEL) { 
     loadstarpar_catalog (stars, Nstars, region, fullname, options);
   } else {
@@ -61 +61 @@
     // got a valid slot, so launch a new host
 
-    // XXX do this after we load the host table:
-    // ensure that the paths are absolute path names
-    // char *tmppath = abspath (table->hosts[i].pathname, DVO_MAX_PATH);
-    // free (table->hosts[i].pathname);
-    // table->hosts[i].pathname = tmppath;
-
     // need to generate the remote command
     char command[1024];