Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/include/relphot.h
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/include/relphot.h	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/include/relphot.h	(revision 31261)
@@ -85,7 +85,9 @@
 int    RELPHOT_GRID_BINNING;
 
-PhotCode      *photcode;
-int            PhotNsec;
-int            PhotSec;
+int      *photseclist;
+int      Nphotcodes;
+PhotCode **photcodes;
+// int            PhotSec;
+// int            PhotNsec;
 
 PhotCode      *refPhotcode;
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/GridOps.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/GridOps.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/GridOps.c	(revision 31261)
@@ -410,5 +410,4 @@
 	// select the color- and airmass-corrected observed magnitude for this star
 	// XXX need to be able to turn off the color-correction until initial average mags are found
-	// Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*PhotNsec]);
 	Msys = PhotCat (&catalog[c].measure[m]);
 	if (isnan(Msys)) {
@@ -528,4 +527,6 @@
   if (!USE_GRID) return;
 
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+
   Nmax = Nlist[0];
   for (i = 0; i < Ngrid; i++) {
@@ -566,5 +567,5 @@
       
       n = catalog[c].measure[m].averef;
-      Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*PhotNsec]);
+      Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*Nsecfilt]);
       if (isnan(Msys)) {
 	Nsys++;
@@ -612,4 +613,6 @@
 
   if (!USE_GRID) return;
+
+  int Nsecfilt = GetPhotcodeNsecfilt ();
 
   N = 0;
@@ -642,5 +645,5 @@
 
       n = catalog[c].measure[m].averef;
-      Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*PhotNsec]);
+      Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*Nsecfilt]);
 
       xlist[N] = Xmeas[c][m];
@@ -704,5 +707,5 @@
   gfits_create_matrix (&header, &matrix);
   gfits_modify (&header, "NEXTEND", OFF_T_FMT, 1,  Nimage + 3);
-  gfits_modify (&header, "FILTER", "%s", 1, photcode[0].name);
+  gfits_modify (&header, "FILTER", "%s", 1, photcodes[0][0].name);  // XXXX note that this expects a single photcode, enforced in initialize.d
   gfits_modify_alt (&header, "COMMENT", "%S", 1, "Mosaic Photometry Grid Analysis");
 
@@ -727,5 +730,5 @@
   theader.bitpix   = -32;
   gfits_create_Theader (&theader, "IMAGE");
-  gfits_modify (&theader, "FILTER", "%s", 1, photcode[0].name);
+  gfits_modify (&theader, "FILTER", "%s", 1, photcodes[0][0].name);
   gfits_modify (&theader, "EXTNAME", "%s", 1, "MAG_OFFSET");
   gfits_create_matrix  (&theader, &matrix);
@@ -742,5 +745,5 @@
   /* save grid Nmeas values */
   gfits_modify (&theader, "EXTNAME", "%s", 1, "NMEAS");
-  gfits_modify (&theader, "FILTER", "%s", 1, photcode[0].name);
+  gfits_modify (&theader, "FILTER", "%s", 1, photcodes[0][0].name);
   gfits_create_matrix  (&theader, &matrix);
   for (i = 0; i < gridX; i++) {
@@ -756,5 +759,5 @@
   /* save grid sigma values */
   gfits_modify (&theader, "EXTNAME", "%s", 1, "SIGMA");
-  gfits_modify (&theader, "FILTER", "%s", 1, photcode[0].name);
+  gfits_modify (&theader, "FILTER", "%s", 1, photcodes[0][0].name);
   gfits_create_matrix  (&theader, &matrix);
   for (i = 0; i < gridX; i++) {
@@ -776,5 +779,5 @@
 
     gfits_modify (&theader, "EXTNAME", "%s", 1, camera.ccdname[N]);
-    gfits_modify (&theader, "FILTER", "%s", 1, photcode[0].name);
+    gfits_modify (&theader, "FILTER", "%s", 1, photcodes[0][0].name);
     gfits_modify (&theader, "NX", "%d", 1, camera.Nx);
     gfits_modify (&theader, "NY", "%d", 1, camera.Ny);
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/ImageOps.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/ImageOps.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/ImageOps.c	(revision 31261)
@@ -128,14 +128,18 @@
 }
 
-/* select all image equivalent to the current photcode */
+/* select all image equivalent to the active photcode set */
 void findImages (Catalog *catalog, int Ncatalog) {
 
   off_t j;
-  int i, ecode;
-
-  for (i = 0; i < Ncatalog; i++) {
+  int i, ecode, Ns, found;
+  
+ for (i = 0; i < Ncatalog; i++) {
     for (j = 0; j < catalog[i].Nmeasure; j++) {
       ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[j].photcode);
-      if (photcode[0].code != ecode) continue;
+      found = FALSE;
+      for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
+	if (ecode == photcodes[Ns][0].code) found = TRUE;
+      }
+      if (!found) continue; 
       matchImage (catalog, j, i);
     }
@@ -145,5 +149,5 @@
 int findCCD (off_t idx, off_t meas, int cat, Measure *measure) {
 
-  int ccdnum;
+  int ccdnum, found, Ns, ecode;
   double X, Y;
   char *pname, *filter, *p, base[256];
@@ -151,5 +155,14 @@
   /* identify the ccd on the basis of the photcode name */
   pname = GetPhotcodeNamebyCode (image[idx].photcode);
-  filter = photcode[0].name;
+
+  // XXX this seems quite terrible...
+  ecode = GetPhotcodeEquivCodebyCode (measure[0].photcode);
+  found = FALSE;
+  for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
+    if (ecode == photcodes[Ns][0].code) found = TRUE;
+  }
+  if (!found) return (FALSE); 
+
+  filter = photcodes[Ns][0].name;
   sprintf (base, "%s.%s.", MOSAICNAME, filter);
   if (strncmp (pname, base, strlen (base))) return (FALSE);
@@ -307,4 +320,6 @@
 
   if (FREEZE_IMAGES) return;
+
+  int Nsecfilt = GetPhotcodeNsecfilt ();
 
   if (PoorImages) {
@@ -360,5 +375,5 @@
       
       n = catalog[c].measure[m].averef;
-      Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*PhotNsec]);
+      Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*Nsecfilt]);
       if (isnan(Msys)) {
 	Nsys++;
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/MosaicOps.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/MosaicOps.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/MosaicOps.c	(revision 31261)
@@ -262,5 +262,5 @@
 int findMosaics (Catalog *catalog, int Ncatalog) {
   
-  int i, ecode;
+  int i, ecode, found, Ns;
   off_t j;
 
@@ -274,5 +274,9 @@
       }
       ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[j].photcode);
-      if (photcode[0].code != ecode) continue;
+      found = FALSE;
+      for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
+	if (ecode == photcodes[Ns][0].code) found = TRUE;
+      }
+      if (!found) continue;
       matchMosaics (catalog, j, i);
     }
@@ -387,5 +391,5 @@
 
   off_t i, j, m, c, N, Nmax;
-  int n, mark, bad, Nfew, Nbad, Ncal, Nrel, Ngrid, Nsys;
+  int n, mark, bad, Nfew, Nbad, Ncal, Nrel, Ngrid, Nsys, Nsecfilt;
   float Msys, Mrel, Mcal, Mgrid;
   double *list, *dlist, *Mlist, *dMlist;
@@ -397,4 +401,6 @@
 
   image = getimages (&N);
+
+  Nsecfilt = GetPhotcodeNsecfilt ();
 
   if (PoorImages) {
@@ -448,5 +454,5 @@
       
       n = catalog[c].measure[m].averef;
-      Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*PhotNsec]);
+      Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*Nsecfilt]);
       if (isnan(Msys)) {
 	Nsys++;
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/StarOps.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/StarOps.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/StarOps.c	(revision 31261)
@@ -29,11 +29,18 @@
 float getMrel (Catalog *catalog, off_t meas, int cat) {
 
+  int Nsec, Nsecfilt, code;
   int ave;
   float value;
 
   ave = catalog[cat].measure[meas].averef;
+  code = catalog[cat].measure[meas].photcode;
+
+  // XXX DEP : this should use a flag associated with secfilt
   if (catalog[cat].average[ave].flags & STAR_BAD) return (NAN);  
  
-  value = catalog[cat].secfilt[PhotNsec*ave+PhotSec].M;
+  Nsec = GetPhotcodeNsec(code);
+  Nsecfilt = GetPhotcodeNsecfilt ();
+
+  value = catalog[cat].secfilt[Nsecfilt*ave+Nsec].M;
   return (value);
 }
@@ -45,5 +52,7 @@
   float Msys, Mcal, Mmos, Mgrid;
   StatType stats;
-
+  int Nsec, Nsecfilt, ecode;
+
+  Nsecfilt = GetPhotcodeNsecfilt ();
   Nfew = Nsys = Nbad = Ncal = Nmos = Ngrid = 0;
 
@@ -51,68 +60,81 @@
     for (j = 0; j < catalog[i].Naverage; j++) {
 
-      /* calculate the average value for a single star */
-      if (catalog[i].average[j].flags & STAR_BAD) continue;  
-      m = catalog[i].average[j].measureOffset;
-
-      N = 0;
-      for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
-	if (catalog[i].measure[m].dbFlags & MEAS_BAD) {
-	  Nbad ++;
-	  continue;
-	}
-	// XXX allow REF stars (no Image Entry) to be included in the calculation this
-	// should be optionally set, and should allow for REF stars to be downweighted by
-	// more than their reported errors.  how such info is carried is unclear...
-	if (getImageEntry (m, i) < 0) {
-	  Mcal = Mmos = Mgrid = 0;
-	} else {
-	  Mcal  = getMcal  (m, i);
-	  if (isnan(Mcal)) { 
-	    Ncal ++;
+      int Ns;
+      for (Ns = 0; Ns < Nphotcodes; Ns++) {
+
+	int thisCode = photcodes[Ns][0].code;
+	Nsec = GetPhotcodeNsec(thisCode);
+
+	/* calculate the average value for a single star */
+	// XXX this flag should be set by secfilt, not average
+	if (catalog[i].average[j].flags & STAR_BAD) continue;  
+	m = catalog[i].average[j].measureOffset;
+
+	N = 0;
+	for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
+	  // skip measurements that do not match the current photcode
+	  ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
+	  if (ecode != thisCode) { continue; }
+
+	  if (catalog[i].measure[m].dbFlags & MEAS_BAD) {
+	    Nbad ++;
 	    continue;
 	  }
-	  Mmos  = getMmos  (m, i);
-	  if (isnan(Mmos)) {
-	    Nmos ++;
+	  // XXX allow REF stars (no Image Entry) to be included in the calculation this
+	  // should be optionally set, and should allow for REF stars to be downweighted by
+	  // more than their reported errors.  how such info is carried is unclear...
+	  if (getImageEntry (m, i) < 0) {
+	    Mcal = Mmos = Mgrid = 0;
+	  } else {
+	    Mcal  = getMcal  (m, i);
+	    if (isnan(Mcal)) { 
+	      Ncal ++;
+	      continue;
+	    }
+	    Mmos  = getMmos  (m, i);
+	    if (isnan(Mmos)) {
+	      Nmos ++;
+	      continue;
+	    }
+	    Mgrid = getMgrid (m, i);
+	    if (isnan(Mgrid)) {
+	      Ngrid++;
+	      continue;
+	    }
+	  }
+
+	  Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
+	  if (isnan(Msys)) {
+	    Nsys++;
 	    continue;
 	  }
-	  Mgrid = getMgrid (m, i);
-	  if (isnan(Mgrid)) {
-	    Ngrid++;
-	    continue;
+	  list[N] = Msys - Mcal - Mmos - Mgrid;
+	  dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
+
+	  // tie down reference photometry if the -refcode (code) option is selected
+	  if (refPhotcode) {
+	    if (GetPhotcodeEquivCodebyCode(catalog[i].measure[m].photcode) == refPhotcode[0].equiv) {
+	      // increase the weight by a factor of 100:
+	      dlist[N] = 0.01*catalog[i].measure[m].dM;
+	    }
 	  }
-	}
-
-	Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*PhotNsec]);
-	if (isnan(Msys)) {
-	  Nsys++;
-	  continue;
-	}
-	list[N] = Msys - Mcal - Mmos - Mgrid;
-	dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
-
-	if (refPhotcode) {
-	  if (GetPhotcodeEquivCodebyCode(catalog[i].measure[m].photcode) == refPhotcode[0].equiv) {
-	    // increase the weight by a factor of 100:
-	    dlist[N] = 0.01*catalog[i].measure[m].dM;
-	  }
-	}
-	N++;
-      }
-
-      // when performing the grid analysis, STAR_TOOFEW will be set to 1;
-
-      if (N <= STAR_TOOFEW) { /* too few measurements */
-	catalog[i].average[j].flags |= ID_STAR_FEW;
-	Nfew ++;
-      } else {
-	catalog[i].average[j].flags &= ~ID_STAR_FEW;
-      }	
-
-      liststats (list, dlist, N, &stats);
-
-      catalog[i].secfilt[PhotNsec*j+PhotSec].M  = stats.mean;
-      catalog[i].secfilt[PhotNsec*j+PhotSec].dM = stats.sigma;
-      catalog[i].secfilt[PhotNsec*j+PhotSec].Xm = (stats.Nmeas > 1) ? 100.0*log10(stats.chisq) : NAN_S_SHORT;
+	  N++;
+	}
+
+	// when performing the grid analysis, STAR_TOOFEW will be set to 1;
+	// XXX DEP this makes no sense: need a separate flag for each secfilt
+	if (N <= STAR_TOOFEW) { /* too few measurements */
+	  catalog[i].average[j].flags |= ID_STAR_FEW;
+	  Nfew ++;
+	} else {
+	  catalog[i].average[j].flags &= ~ID_STAR_FEW;
+	}	
+
+	liststats (list, dlist, N, &stats);
+
+	catalog[i].secfilt[Nsecfilt*j+Nsec].M  = stats.mean;
+	catalog[i].secfilt[Nsecfilt*j+Nsec].dM = stats.sigma;
+	catalog[i].secfilt[Nsecfilt*j+Nsec].Xm = (stats.Nmeas > 1) ? 100.0*log10(stats.chisq) : NAN_S_SHORT;
+      }
     }
   }
@@ -129,4 +151,7 @@
   double *list, *dlist;
   StatType stats;
+  int Nsec, Nsecfilt, ecode;
+
+  Nsecfilt = GetPhotcodeNsecfilt ();
 
   /* Nmeasure is now different, need to reallocate */
@@ -142,40 +167,50 @@
   for (i = 0; i < Ncatalog; i++) {
     for (j = 0; j < catalog[i].Naverage; j++) {
-
       /* skip stars already calibrated */
       if (catalog[i].found[j]) continue;  
 
-      N = 0;
-      m = catalog[i].average[j].measureOffset;
-      for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
-	if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue;
-	// XXX allow REF stars (no Image Entry) to be included in the calculation this
-	// should be optionally set, and should allow for REF stars to be downweighted by
-	// more than their reported errors.  how such info is carried is unclear...
-	if (getImageEntry (m, i) < 0) {
-	  Mcal = Mmos = Mgrid = 0;
-	} else {
-	  Mcal  = getMcal  (m, i);
-	  if (isnan(Mcal)) continue;
-	  Mmos  = getMmos  (m, i);
-	  if (isnan(Mmos)) continue;
-	  Mgrid = getMgrid (m, i);
-	  if (isnan(Mgrid)) continue;
-	}
-
-	Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*PhotNsec]);
-	list[N] = Msys - Mcal - Mmos - Mgrid;
-	dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
-	N++;
-      }
-      if (N < 1) continue;
-
-      liststats (list, dlist, N, &stats);
-      if (mark) catalog[i].found[j] = TRUE;
-
-      /* use sigma or error in dM for output? */
-      catalog[i].secfilt[PhotNsec*j+PhotSec].M  = stats.mean;
-      catalog[i].secfilt[PhotNsec*j+PhotSec].dM = MAX (stats.error, stats.sigma);
-      catalog[i].secfilt[PhotNsec*j+PhotSec].Xm = (stats.Nmeas > 1) ? 100.0*log10(stats.chisq) : NAN_S_SHORT;
+      int Ns;
+      for (Ns = 0; Ns < Nphotcodes; Ns++) {
+	int thisCode = photcodes[Ns][0].code;
+	Nsec = GetPhotcodeNsec(thisCode);
+
+	N = 0;
+	m = catalog[i].average[j].measureOffset;
+	for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
+	  // skip measurements that do not match the current photcode
+	  ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
+	  if (ecode != thisCode) { continue; }
+
+	  if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue;
+
+	  // XXX allow REF stars (no Image Entry) to be included in the calculation this
+	  // should be optionally set, and should allow for REF stars to be downweighted by
+	  // more than their reported errors.  how such info is carried is unclear...
+	  if (getImageEntry (m, i) < 0) {
+	    Mcal = Mmos = Mgrid = 0;
+	  } else {
+	    Mcal  = getMcal  (m, i);
+	    if (isnan(Mcal)) continue;
+	    Mmos  = getMmos  (m, i);
+	    if (isnan(Mmos)) continue;
+	    Mgrid = getMgrid (m, i);
+	    if (isnan(Mgrid)) continue;
+	  }
+
+	  Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
+	  list[N] = Msys - Mcal - Mmos - Mgrid;
+	  dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
+	  N++;
+	}
+	if (N < 1) continue;
+
+	liststats (list, dlist, N, &stats);
+	if (mark) catalog[i].found[j] = TRUE;
+
+	/* use sigma or error in dM for output? */
+	catalog[i].secfilt[Nsecfilt*j+Nsec].M  = stats.mean;
+	catalog[i].secfilt[Nsecfilt*j+Nsec].dM = MAX (stats.error, stats.sigma);
+	catalog[i].secfilt[Nsecfilt*j+Nsec].Xm = (stats.Nmeas > 1) ? 100.0*log10(stats.chisq) : NAN_S_SHORT;
+      }
     }
   }
@@ -209,5 +244,5 @@
   ALLOCATE (dlist, double, MAX (1, Nmax));
 
-  Nsecfilt = catalog[0].Nsecfilt;
+  Nsecfilt = GetPhotcodeNsecfilt ();
 
 # define PSFQUALSTATS 1
@@ -227,5 +262,5 @@
 	  if (GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode) != Nc) continue;
 
-	  Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*PhotNsec]);
+	  Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
 	  if (isnan(Msys)) continue;
 
@@ -346,5 +381,5 @@
 void clean_stars (Catalog *catalog, int Ncatalog) {
 
-  int i, j, Ndel, Nave, Ntot, mark;
+  int i, j, Ndel, Nave, Ntot, mark, Ns;
   float dM, Xm;
   double Chisq, MaxScatter, MaxChisq;
@@ -361,42 +396,53 @@
   ALLOCATE (slist, double, Ntot);
   ALLOCATE (dlist, double, Ntot);
-  for (i = Ntot = 0; i < Ncatalog; i++) {
-    for (j = 0; j < catalog[i].Naverage; j++) {
-      if (catalog[i].average[j].flags & STAR_BAD) continue;
-      Xm = catalog[i].secfilt[PhotNsec*j+PhotSec].Xm;
-      if (Xm == -1) continue;
-      Chisq = pow (10.0, 0.01*Xm);
-      xlist[Ntot] = Chisq;
-      slist[Ntot] = catalog[i].secfilt[PhotNsec*j+PhotSec].dM;
-      dlist[Ntot] = 1;
-      Ntot ++;
-    }
-  }
+
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+
+  // eliminate bad stars using the stats for a single secfilt at a time
+  // XXX DEP replace average.flags with secfilt flags
+  for (Ns = 0; Ns < Nphotcodes; Ns ++) {
+    
+    int thisCode = photcodes[Ns][0].code;
+    int Nsec = GetPhotcodeNsec(thisCode);
+
+    for (i = Ntot = 0; i < Ncatalog; i++) {
+      for (j = 0; j < catalog[i].Naverage; j++) {
+	if (catalog[i].average[j].flags & STAR_BAD) continue;
+	Xm = catalog[i].secfilt[Nsecfilt*j+Nsec].Xm;
+	if (Xm == -1) continue;
+	Chisq = pow (10.0, 0.01*Xm);
+	xlist[Ntot] = Chisq;
+	slist[Ntot] = catalog[i].secfilt[Nsecfilt*j+Nsec].dM;
+	dlist[Ntot] = 1;
+	Ntot ++;
+      }
+    }
   
-  initstats ("MEAN");
-  liststats (xlist, dlist, Ntot, &stats);
-  MaxChisq = MAX (STAR_CHISQ, 2*stats.median);
-  liststats (slist, dlist, Ntot, &stats);
-  MaxScatter = MAX (STAR_SCATTER, 2*stats.median);
-  fprintf (stderr, "Max Scatter: %f, Max Chisq: %f\n", MaxScatter, MaxChisq);
-
-  Ndel = Nave = 0;
-  for (i = 0; i < Ncatalog; i++) {
-    for (j = 0; j < catalog[i].Naverage; j++) {
-      dM = catalog[i].secfilt[PhotNsec*j+PhotSec].dM;
-      Xm = catalog[i].secfilt[PhotNsec*j+PhotSec].Xm;
-      Chisq = pow (10.0, 0.01*Xm);
-      mark = (dM > MaxScatter) || (Xm == NAN_S_SHORT) || (Chisq > MaxChisq);
-      if (mark) {
-	catalog[i].average[j].flags |= ID_STAR_POOR;
-	Ndel ++;
-      } else {
-	catalog[i].average[j].flags &= ~ID_STAR_POOR;
-      }
-      Nave ++;
-    }
-  }
-  fprintf (stderr, "%d stars marked variable, %d total\n", Ndel, Nave);
-  initstats (STATMODE);
+    initstats ("MEAN");
+    liststats (xlist, dlist, Ntot, &stats);
+    MaxChisq = MAX (STAR_CHISQ, 2*stats.median);
+    liststats (slist, dlist, Ntot, &stats);
+    MaxScatter = MAX (STAR_SCATTER, 2*stats.median);
+    fprintf (stderr, "Max Scatter: %f, Max Chisq: %f\n", MaxScatter, MaxChisq);
+
+    Ndel = Nave = 0;
+    for (i = 0; i < Ncatalog; i++) {
+      for (j = 0; j < catalog[i].Naverage; j++) {
+	dM = catalog[i].secfilt[Nsecfilt*j+Nsec].dM;
+	Xm = catalog[i].secfilt[Nsecfilt*j+Nsec].Xm;
+	Chisq = pow (10.0, 0.01*Xm);
+	mark = (dM > MaxScatter) || (Xm == NAN_S_SHORT) || (Chisq > MaxChisq);
+	if (mark) {
+	  catalog[i].average[j].flags |= ID_STAR_POOR;
+	  Ndel ++;
+	} else {
+	  catalog[i].average[j].flags &= ~ID_STAR_POOR;
+	}
+	Nave ++;
+      }
+    }
+    fprintf (stderr, "%d stars marked variable, %d total\n", Ndel, Nave);
+    initstats (STATMODE);
+  }
   free (xlist);
   free (slist);
@@ -416,4 +462,6 @@
   int Ncal, Nmos, Ngrid, Nfew;
 
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+
   if (VERBOSE) fprintf (stderr, "marking poor measures\n");
   /* Nmeasure is now different, need to reallocate */
@@ -437,79 +485,93 @@
     for (j = 0; j < catalog[i].Naverage; j++) {
 
-      /* skip bad stars to prevent them from becoming good (on inner sample) */
-      if (catalog[i].average[j].flags & STAR_BAD) continue;  
-
-      /* on final processing, skip stars already measured */
-      if (final && catalog[i].found[j]) continue;  
-
-      /* accumulate list of valid measurements */
-      m = catalog[i].average[j].measureOffset;
-      N = 0;
-      for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
-	/* if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue; */
-	Mcal  = getMcal  (m, i);
-	if (isnan(Mcal)) { Ncal ++; continue; }
-	Mmos  = getMmos  (m, i);
-	if (isnan(Mmos)) { Nmos ++; continue; }
-	Mgrid = getMgrid (m, i);
-	if (isnan(Mgrid)) { Ngrid ++; continue; }
-
-	Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*PhotNsec]);
-	list[N] = Msys - Mcal - Mmos - Mgrid;
-	dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
-	N++;
-      }
-      if (N <= TOOFEW) { Nfew ++; continue; }
-
-      /* 3-sigma clip based on stats of inner 50% */
-
-      // calculated mean of inner 50%
-      initstats ("INNER_MEAN");
-      liststats (list, dlist, N, &stats);
-      stats.sigma = MAX (MIN_ERROR, stats.sigma); /* if measurements agree too well, sigma -> 0.0 */
-
-      // ignore entries > 3sigma from inner mean
-      for (k = m = 0; k < N; k++) {
-	if (fabs (list[k] - stats.median) < NSIGMA_CLIP*stats.sigma) {
-	  list[m] = list[k];
-	  m++;
-	}
-      }
-      // recalculate the mean & sigma of the accepted measurements
-      initstats ("MEAN");
-      liststats (list, dlist, m, &stats);
-      stats.sigma = MAX (MIN_ERROR, stats.sigma);
-
-      /* apply to list of all relevant measurements, including IMAGE_POOR & IMAGE_FEW */
-      image_bad = IMAGE_BAD;
-      IMAGE_BAD = ID_IMAGE_PHOTOM_NOCAL;
-      m = catalog[i].average[j].measureOffset;
-      N = 0;
-      for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
-	/* if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue; */
-	Mcal  = getMcal  (m, i);
-	if (isnan(Mcal)) continue;
-	Mmos  = getMmos  (m, i);
-	if (isnan(Mmos)) continue;
-	Mgrid = getMgrid (m, i);
-	if (isnan(Mgrid)) continue;
-
-	Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*PhotNsec]);
-	list[N] = Msys - Mcal - Mmos - Mgrid;
-	dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
-	ilist[N] = m;
-	N++;
-	Nave ++;
-      }
-      if (N < TOOFEW) continue;
-
-      /* mark bad measures (> 3 sigma deviant) */
-      for (k = 0; k < N; k++) {
-	if (fabs (list[k] - stats.median) > NSIGMA_REJECT*stats.sigma) {
-	  catalog[i].measure[ilist[k]].dbFlags |= ID_MEAS_POOR_PHOTOM;
-	  Ndel ++;
-	}
-      }
-      IMAGE_BAD = image_bad;
+      int Ns;
+      for (Ns = 0; Ns < Nphotcodes; Ns++) {
+	
+	/* on final processing, skip stars already measured */
+	if (final && catalog[i].found[j]) continue;  
+
+	int thisCode = photcodes[Ns][0].code;
+
+	/* skip bad stars to prevent them from becoming good (on inner sample) */
+	if (catalog[i].average[j].flags & STAR_BAD) continue;  
+
+	/* accumulate list of valid measurements */
+	m = catalog[i].average[j].measureOffset;
+	N = 0;
+	for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
+	  // skip measurements that do not match the current photcode
+	  int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
+	  if (ecode != thisCode) { continue; }
+
+	  /* if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue; */
+	  Mcal  = getMcal  (m, i);
+	  if (isnan(Mcal)) { Ncal ++; continue; }
+	  Mmos  = getMmos  (m, i);
+	  if (isnan(Mmos)) { Nmos ++; continue; }
+	  Mgrid = getMgrid (m, i);
+	  if (isnan(Mgrid)) { Ngrid ++; continue; }
+
+	  Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
+	  list[N] = Msys - Mcal - Mmos - Mgrid;
+	  dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
+	  N++;
+	}
+	if (N <= TOOFEW) { Nfew ++; continue; }
+
+	/* 3-sigma clip based on stats of inner 50% */
+
+	// calculated mean of inner 50%
+	initstats ("INNER_MEAN");
+	liststats (list, dlist, N, &stats);
+	stats.sigma = MAX (MIN_ERROR, stats.sigma); /* if measurements agree too well, sigma -> 0.0 */
+
+	// ignore entries > 3sigma from inner mean
+	for (k = m = 0; k < N; k++) {
+	  if (fabs (list[k] - stats.median) < NSIGMA_CLIP*stats.sigma) {
+	    list[m] = list[k];
+	    m++;
+	  }
+	}
+	// recalculate the mean & sigma of the accepted measurements
+	initstats ("MEAN");
+	liststats (list, dlist, m, &stats);
+	stats.sigma = MAX (MIN_ERROR, stats.sigma);
+
+	/* apply to list of all relevant measurements, including IMAGE_POOR & IMAGE_FEW */
+	image_bad = IMAGE_BAD;
+	IMAGE_BAD = ID_IMAGE_PHOTOM_NOCAL;
+	m = catalog[i].average[j].measureOffset;
+	N = 0;
+	for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
+	  // skip measurements that do not match the current photcode
+	  int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
+	  if (ecode != thisCode) { continue; }
+
+	  /* if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue; */
+	  Mcal  = getMcal  (m, i);
+	  if (isnan(Mcal)) continue;
+	  Mmos  = getMmos  (m, i);
+	  if (isnan(Mmos)) continue;
+	  Mgrid = getMgrid (m, i);
+	  if (isnan(Mgrid)) continue;
+
+	  Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
+	  list[N] = Msys - Mcal - Mmos - Mgrid;
+	  dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
+	  ilist[N] = m;
+	  N++;
+	  Nave ++;
+	}
+	if (N < TOOFEW) continue;
+
+	/* mark bad measures (> 3 sigma deviant) */
+	for (k = 0; k < N; k++) {
+	  if (fabs (list[k] - stats.median) > NSIGMA_REJECT*stats.sigma) {
+	    catalog[i].measure[ilist[k]].dbFlags |= ID_MEAS_POOR_PHOTOM;
+	    Ndel ++;
+	  }
+	}
+	IMAGE_BAD = image_bad;
+      }
     }
   }
@@ -567,4 +629,5 @@
 }
 
+// XXXX do this for all secfilt at once or one at a time?
 StatType statsStarX (Catalog *catalog, int Ncatalog) {
 
@@ -574,4 +637,6 @@
   StatType stats;
 
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+
   Ntot = 0;
   for (i = 0; i < Ncatalog; i++) {
@@ -586,12 +651,19 @@
     for (j = 0; j < catalog[i].Naverage; j++) {
 
-      /* calculate the average value for a single star */
-      if (catalog[i].average[j].flags & STAR_BAD) continue;  
-
-      Xm = catalog[i].secfilt[PhotNsec*j+PhotSec].Xm;
-      if (Xm == NAN_S_SHORT) continue;
-      list[n] = pow (10.0, 0.01*Xm);
-      dlist[n] = 1;
-      n++;
+      int Ns;
+      for (Ns = 0; Ns < Nphotcodes; Ns++) {
+
+	int thisCode = photcodes[Ns][0].code;
+	int Nsec = GetPhotcodeNsec(thisCode);
+
+	/* calculate the average value for a single star */
+	if (catalog[i].average[j].flags & STAR_BAD) continue;  
+	
+	Xm = catalog[i].secfilt[Nsecfilt*j+Nsec].Xm;
+	if (Xm == NAN_S_SHORT) continue;
+	list[n] = pow (10.0, 0.01*Xm);
+	dlist[n] = 1;
+	n++;
+      }
     }
   }
@@ -611,4 +683,6 @@
   StatType stats;
 
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+
   Ntot = 0;
   for (i = 0; i < Ncatalog; i++) {
@@ -623,11 +697,18 @@
     for (j = 0; j < catalog[i].Naverage; j++) {
 
-      /* calculate the average value for a single star */
-      if (catalog[i].average[j].flags & STAR_BAD) continue;  
-
-      dM = catalog[i].secfilt[PhotNsec*j+PhotSec].dM;
-      list[n] = dM;
-      dlist[n] = 1;
-      n++;
+      int Ns;
+      for (Ns = 0; Ns < Nphotcodes; Ns++) {
+
+	int thisCode = photcodes[Ns][0].code;
+	int Nsec = GetPhotcodeNsec(thisCode);
+
+	/* calculate the average value for a single star */
+	if (catalog[i].average[j].flags & STAR_BAD) continue;  
+
+	dM = catalog[i].secfilt[Nsecfilt*j+Nsec].dM;
+	list[n] = dM;
+	dlist[n] = 1;
+	n++;
+      }
     }
   }
@@ -651,20 +732,28 @@
   ALLOCATE (Mlist, double, NBIN);
 
-  for (i = 0; i < NBIN; i++) xlist[i] = 0.00025*i;
-  bzero (Mlist, NBIN*sizeof(double));
-  for (i = 0; i < Ncatalog; i++) {
-    for (j = 0; j < catalog[i].Naverage; j++) {
-      if (catalog[i].average[j].flags & STAR_BAD) continue;  
-      dMrel = catalog[i].secfilt[PhotNsec*j+PhotSec].dM;
-      bin = dMrel / 0.00025;
-      bin = MAX (0, MIN (NBIN-1, bin));
-      Mlist[bin] += 1.0;
-    }
-  }
-
-  plot_defaults (&graphdata);
-  graphdata.style = 1;
-  plot_list (&graphdata, xlist, Mlist, NBIN, "dMrel hist", "%s.dMhist.png", OUTROOT);
-
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+
+  int Ns;
+  for (Ns = 0; Ns < Nphotcodes; Ns++) {
+
+    int thisCode = photcodes[Ns][0].code;
+    int Nsec = GetPhotcodeNsec(thisCode);
+
+    for (i = 0; i < NBIN; i++) xlist[i] = 0.00025*i;
+    bzero (Mlist, NBIN*sizeof(double));
+    for (i = 0; i < Ncatalog; i++) {
+      for (j = 0; j < catalog[i].Naverage; j++) {
+	if (catalog[i].average[j].flags & STAR_BAD) continue;  
+	dMrel = catalog[i].secfilt[Nsecfilt*j+Nsec].dM;
+	bin = dMrel / 0.00025;
+	bin = MAX (0, MIN (NBIN-1, bin));
+	Mlist[bin] += 1.0;
+      }
+    }
+
+    plot_defaults (&graphdata);
+    graphdata.style = 1;
+    plot_list (&graphdata, xlist, Mlist, NBIN, "dMrel hist", "%s.dMhist.png", OUTROOT);
+  }
   free (xlist);
   free (Mlist);
@@ -678,4 +767,6 @@
   Graphdata graphdata;
 
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+
   Ntotal = 0;
   for (i = 0; i < Ncatalog; i++) Ntotal += catalog[i].Naverage;
@@ -684,19 +775,27 @@
   ALLOCATE (ylist, double, Ntotal);
 
-  N = 0;
-  for (i = 0; i < Ncatalog; i++) {
-    for (j = 0; j < catalog[i].Naverage; j++) {
-      if (catalog[i].average[j].flags & STAR_BAD) continue;
-      xlist[N] = catalog[i].secfilt[PhotNsec*j+PhotSec].M;
-      value    = catalog[i].secfilt[PhotNsec*j+PhotSec].Xm;
-      if (value == NAN_S_SHORT) continue;
-      ylist[N] = 0.01*value;
-      N++;
-    }
-  }
-
-  plot_defaults (&graphdata);
-  graphdata.ymin = -3.0;
-  plot_list (&graphdata, xlist, ylist, N, "chisq", "%s.chisq.png", OUTROOT);
+  int Ns;
+  for (Ns = 0; Ns < Nphotcodes; Ns++) {
+
+    int thisCode = photcodes[Ns][0].code;
+    int Nsec = GetPhotcodeNsec(thisCode);
+
+    N = 0;
+    for (i = 0; i < Ncatalog; i++) {
+      for (j = 0; j < catalog[i].Naverage; j++) {
+	if (catalog[i].average[j].flags & STAR_BAD) continue;
+	xlist[N] = catalog[i].secfilt[Nsecfilt*j+Nsec].M;
+	value    = catalog[i].secfilt[Nsecfilt*j+Nsec].Xm;
+	if (value == NAN_S_SHORT) continue;
+	ylist[N] = 0.01*value;
+	N++;
+      }
+    }
+
+    plot_defaults (&graphdata);
+    graphdata.ymin = -3.0;
+    plot_list (&graphdata, xlist, ylist, N, "chisq", "%s.chisq.png", OUTROOT);
+  }
+
   free (xlist);
   free (ylist);
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/bcatalog.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/bcatalog.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/bcatalog.c	(revision 31261)
@@ -6,13 +6,11 @@
   
   off_t i, j, offset;
-  int ecode;
+  int ecode, found, Ns;
   off_t NAVERAGE, NMEASURE, Naverage, Nmeasure, Nm;
   float mag;
   int Ncode, Ntime, Ndophot, Nmag, Nsigma, Nimag, Nfew, Ngalaxy, Npsfqf;
 
-  // XXX PhotNsec as a global is a bad idea; either get it from catalog
-  // or get it from:
-  // Nsecfilt = GetPhotcodeNsecfilt ();
-  // assert (catalog[0].Nsecfilt == Nsecfilt);
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+  assert (Nsecfilt == catalog[0].Nsecfilt);
 
   /* we are moving only the subset of measurements from catalog[0] to subcatalog[0] */
@@ -21,5 +19,5 @@
   NMEASURE = 1000;
   ALLOCATE (subcatalog[0].average, Average, NAVERAGE);
-  ALLOCATE (subcatalog[0].secfilt, SecFilt, NAVERAGE*PhotNsec);
+  ALLOCATE (subcatalog[0].secfilt, SecFilt, NAVERAGE*Nsecfilt);
   ALLOCATE (subcatalog[0].measure, Measure, NMEASURE);
   Nmeasure = Naverage = 0;
@@ -34,13 +32,20 @@
     subcatalog[0].average[Naverage] = catalog[0].average[i];
     subcatalog[0].average[Naverage].measureOffset = Nmeasure;
-    for (j = 0; j < PhotNsec; j++) {
-      subcatalog[0].secfilt[PhotNsec*Naverage+j] = catalog[0].secfilt[PhotNsec*i+j];
+    for (j = 0; j < Nsecfilt; j++) {
+      subcatalog[0].secfilt[Nsecfilt*Naverage+j] = catalog[0].secfilt[Nsecfilt*i+j];
     }
 
     if (RESET) {
-      subcatalog[0].secfilt[PhotNsec*Naverage+PhotSec].M  = NAN;
-      subcatalog[0].secfilt[PhotNsec*Naverage+PhotSec].dM = NAN;
-      subcatalog[0].average[Naverage].flags &= ~ID_STAR_FEW;
-      subcatalog[0].average[Naverage].flags &= ~ID_STAR_POOR;
+      int Ns;
+      for (Ns = 0; Ns < Nphotcodes; Ns++) {
+
+	int thisCode = photcodes[Ns][0].code;
+	int Nsec = GetPhotcodeNsec(thisCode);
+
+	subcatalog[0].secfilt[Nsecfilt*Naverage+Nsec].M  = NAN;
+	subcatalog[0].secfilt[Nsecfilt*Naverage+Nsec].dM = NAN;
+	subcatalog[0].average[Naverage].flags &= ~ID_STAR_FEW;
+	subcatalog[0].average[Naverage].flags &= ~ID_STAR_POOR;
+      }
     }
 
@@ -54,5 +59,12 @@
       /* select measurements by photcode */
       ecode = GetPhotcodeEquivCodebyCode (catalog[0].measure[offset].photcode);
-      if (ecode != photcode[0].code) { Ncode ++; continue; }
+      found = FALSE;
+      for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
+	if (ecode == photcodes[Ns][0].code) found = TRUE;
+      }
+      if (!found) {
+	Ncode ++; 
+	continue; 
+      }
 
       /* select measurements by time */
@@ -129,15 +141,14 @@
       NAVERAGE += 50;
       REALLOCATE (subcatalog[0].average, Average, NAVERAGE);
-      REALLOCATE (subcatalog[0].secfilt, SecFilt, NAVERAGE*PhotNsec);
+      REALLOCATE (subcatalog[0].secfilt, SecFilt, NAVERAGE*Nsecfilt);
     }
   }
   REALLOCATE (subcatalog[0].average, Average, MAX (Naverage, 1));
   REALLOCATE (subcatalog[0].measure, Measure, MAX (Nmeasure, 1));
-  REALLOCATE (subcatalog[0].secfilt, SecFilt, PhotNsec*MAX (Naverage, 1));
+  REALLOCATE (subcatalog[0].secfilt, SecFilt, Nsecfilt*MAX (Naverage, 1));
   subcatalog[0].Naverage = Naverage;
   subcatalog[0].Nmeasure = Nmeasure;
   subcatalog[0].Nsecfilt = catalog[0].Nsecfilt;
   subcatalog[0].Nsecf_mem = Naverage * catalog[0].Nsecfilt;
-  assert (PhotNsec == catalog[0].Nsecfilt);
 
   // limit the total number of stars in the catalog
@@ -161,4 +172,6 @@
   double Rmin, Rmax, Dmin, Dmax;
 
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+
   gfits_scan (&catalog[0].header, "RA0",  "%lf", 1, &Rmin);
   gfits_scan (&catalog[0].header, "DEC0", "%lf", 1, &Dmin);
@@ -205,5 +218,5 @@
   // allocate the output data 
   ALLOCATE (tmpcatalog.average, Average, Nmax);
-  ALLOCATE (tmpcatalog.secfilt, SecFilt, Nmax * PhotNsec);
+  ALLOCATE (tmpcatalog.secfilt, SecFilt, Nmax * Nsecfilt);
   ALLOCATE (tmpcatalog.measure, Measure, NMEASURE);
 
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/initialize.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/initialize.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/initialize.c	(revision 31261)
@@ -13,4 +13,6 @@
     N = UserPatchSelect ? 1 : 2;
 
+# if (0)
+    // XXX DEP
     if ((photcode = GetPhotcodebyName (argv[N])) == NULL) {
       fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", argv[N]);
@@ -21,8 +23,39 @@
       exit (1);
     }
+    // PhotSec is used to select the single average photcode being processed
     PhotSec = GetPhotcodeNsec (photcode[0].code);
+# endif
+
+    Nphotcodes = 0;
+    photcodes = NULL;
+    int NPHOTCODES = 10;
+    ALLOCATE (photcodes, PhotCode *, NPHOTCODES);
+
+    /* parse the comma-separated list of photcodesKeep */
+    char *myList = strcreate(argv[N]);
+    char *list = myList;
+    char *codename = NULL;
+    char *ptr = NULL;
+    while ((codename = strtok_r (list, ",", &ptr)) != NULL) {
+      list = NULL; // pass NULL on successive strtok_r calls
+      fprintf (stderr, "PHOTCODE LIST: %s\n", myList);
+      fprintf (stderr, "codename: %s\n", codename);
+      if ((photcodes[Nphotcodes] = GetPhotcodebyName (codename)) == NULL) {
+        fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", codename);
+        exit (1);
+      }
+      if (photcodes[Nphotcodes][0].type != PHOT_SEC) {
+	  fprintf (stderr, "photcode %s is not an filter type (SEC)\n", codename);
+	  exit (1);
+      }
+      Nphotcodes ++;
+      CHECK_REALLOCATE (photcodes, PhotCode *, NPHOTCODES, Nphotcodes, 10);
+    }
   }
-
-  PhotNsec = GetPhotcodeNsecfilt ();
+  // XXX DEP PhotNsec = GetPhotcodeNsecfilt ();
+  if (USE_GRID && (Nphotcodes > 1)) {
+    fprintf (stderr, "grid correction analysis currently can only operate on a single photcode\n");
+    exit (1);
+  }
 
   initstats (STATMODE);
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/plot_scatter.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/plot_scatter.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/plot_scatter.c	(revision 31261)
@@ -18,41 +18,54 @@
   ALLOCATE (ilist, double, Ntot);
 
-  N = 0;
-  for (i = 0; i < Ncatalog; i++) {
-    for (j = 0; j < catalog[i].Naverage; j++) {
+  int Nsecfilt = GetPhotcodeNsecfilt ();
 
-      /* calculate the average value for a single star */
-      if (catalog[i].average[j].flags & STAR_BAD) continue;  
-      m = catalog[i].average[j].measureOffset;
+  int Ns;
+  for (Ns = 0; Ns < Nphotcodes; Ns++) {
 
-      for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
-	if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue;
-	Mcal = getMcal  (m, i);
-	if (isnan(Mcal)) continue;
-	Mmos = getMmos  (m, i);
-	if (isnan(Mmos)) continue;
-	Mgrid = getMgrid (m, i);
-	if (isnan(Mgrid)) continue;
+    int thisCode = photcodes[Ns][0].code;
+    int Nsec = GetPhotcodeNsec(thisCode);
 
-	Mrel = catalog[i].secfilt[PhotNsec*j+PhotSec].M;
-	xlist[N] = Mrel;
-	ylist[N] = PhotSys  (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*PhotNsec]) - Mcal - Mmos - Mgrid - Mrel;
-	ilist[N] = PhotInst (&catalog[i].measure[m]);
-	N++;
-      }
+    N = 0;
+    for (i = 0; i < Ncatalog; i++) {
+	for (j = 0; j < catalog[i].Naverage; j++) {
+
+	    /* calculate the average value for a single star */
+	    if (catalog[i].average[j].flags & STAR_BAD) continue;  
+	    m = catalog[i].average[j].measureOffset;
+
+	    for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
+		// skip measurements that do not match the current photcode
+		int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
+		if (ecode != thisCode) { continue; }
+
+		if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue;
+		Mcal = getMcal  (m, i);
+		if (isnan(Mcal)) continue;
+		Mmos = getMmos  (m, i);
+		if (isnan(Mmos)) continue;
+		Mgrid = getMgrid (m, i);
+		if (isnan(Mgrid)) continue;
+
+		Mrel = catalog[i].secfilt[Nsecfilt*j+Nsec].M;
+		xlist[N] = Mrel;
+		ylist[N] = PhotSys  (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]) - Mcal - Mmos - Mgrid - Mrel;
+		ilist[N] = PhotInst (&catalog[i].measure[m]);
+		N++;
+	    }
+	}
     }
+
+    plot_defaults (&graphdata);
+    graphdata.xmin = PlotMmin;
+    graphdata.xmax = PlotMmax;
+    graphdata.ymin = PlotdMmin;
+    graphdata.ymax = PlotdMmax;
+    plot_list (&graphdata, xlist, ylist, N, "mag vs dmag", "%s.Mag.png", OUTROOT);
+
+    plot_defaults (&graphdata);
+    graphdata.ymin = PlotdMmin;
+    graphdata.ymax = PlotdMmax;
+    plot_list (&graphdata, ilist, ylist, N, "imag vs dmag", "%s.iMag.png", OUTROOT);
   }
-
-  plot_defaults (&graphdata);
-  graphdata.xmin = PlotMmin;
-  graphdata.xmax = PlotMmax;
-  graphdata.ymin = PlotdMmin;
-  graphdata.ymax = PlotdMmax;
-  plot_list (&graphdata, xlist, ylist, N, "mag vs dmag", "%s.Mag.png", OUTROOT);
-
-  plot_defaults (&graphdata);
-  graphdata.ymin = PlotdMmin;
-  graphdata.ymax = PlotdMmax;
-  plot_list (&graphdata, ilist, ylist, N, "imag vs dmag", "%s.iMag.png", OUTROOT);
   free (xlist);
   free (ylist);
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/select_images.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/select_images.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/select_images.c	(revision 31261)
@@ -119,5 +119,10 @@
     /* exclude images by photcode */
     ecode = GetPhotcodeEquivCodebyCode (timage[i].photcode);
-    if (ecode != photcode[0].code) continue;
+    found = FALSE;
+    int Ns;
+    for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
+	if (ecode == photcodes[Ns][0].code) found = TRUE;
+    }
+    if (!found) continue;
 
     /* exclude images by time */
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/setExclusions.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/setExclusions.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/setExclusions.c	(revision 31261)
@@ -4,5 +4,5 @@
 
   off_t i, j, k, m, Narea, Nnocal, Ngood;
-  int ecode;
+  int ecode, found, Ns;
   Coords *coords;
   double r, d, x, y;
@@ -16,5 +16,9 @@
 	/* select measurements by photcode */
 	ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
-	if (ecode != photcode[0].code) goto mark_nocal;
+	found = FALSE;
+	for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
+	    if (ecode == photcodes[Ns][0].code) found = TRUE;
+	}
+	if (!found) goto mark_nocal;
 	
 	/* select measurements by time */
Index: branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/setMrelFinal.c
===================================================================
--- branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/setMrelFinal.c	(revision 31260)
+++ branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/setMrelFinal.c	(revision 31261)
@@ -9,28 +9,37 @@
   if (RESET) {
 
-    for (i = 0; i < catalog[0].Naverage; i++) {
-      catalog[0].secfilt[PhotNsec*i+PhotSec].M  = NAN;
-      catalog[0].secfilt[PhotNsec*i+PhotSec].dM = NAN;
-      catalog[0].secfilt[PhotNsec*i+PhotSec].Xm = NAN_S_SHORT;
+    int Nsecfilt = GetPhotcodeNsecfilt ();
 
-      m = catalog[0].average[i].measureOffset;
-      for (j = 0; j < catalog[0].average[i].Nmeasure; j++, m++) {
+    int Ns;
+    for (Ns = 0; Ns < Nphotcodes; Ns++) {
 	
-	/* select measurements by photcode */
-	ecode = GetPhotcodeEquivCodebyCode (catalog[0].measure[m].photcode);
-	if (ecode != photcode[0].code) continue;
+      int thisCode = photcodes[Ns][0].code;
+      int Nsec = GetPhotcodeNsec(thisCode);
+
+      for (i = 0; i < catalog[0].Naverage; i++) {
+	catalog[0].secfilt[Nsecfilt*i+Nsec].M  = NAN;
+	catalog[0].secfilt[Nsecfilt*i+Nsec].dM = NAN;
+	catalog[0].secfilt[Nsecfilt*i+Nsec].Xm = NAN_S_SHORT;
+
+	m = catalog[0].average[i].measureOffset;
+	for (j = 0; j < catalog[0].average[i].Nmeasure; j++, m++) {
 	
-	/* select measurements by time */
-	if (TimeSelect) {
-	  if (catalog[0].measure[m].t < TSTART) continue;
-	  if (catalog[0].measure[m].t > TSTOP) continue;
+	  // skip measurements that do not match the current photcode
+	  ecode = GetPhotcodeEquivCodebyCode (catalog[0].measure[m].photcode);
+	  if (ecode != thisCode) { continue; }
+
+	  /* select measurements by time */
+	  if (TimeSelect) {
+	    if (catalog[0].measure[m].t < TSTART) continue;
+	    if (catalog[0].measure[m].t > TSTOP) continue;
+	  }
+	
+	  catalog[0].measure[m].Mcal = 0;
+	  catalog[0].measure[m].dbFlags &= 0xff00;
+	  catalog[0].measure[m].dbFlags &= ~ID_MEAS_POOR_PHOTOM;
+	  catalog[0].measure[m].dbFlags &= ~ID_MEAS_SKIP_PHOTOM;
+	  catalog[0].measure[m].dbFlags &= ~ID_MEAS_AREA;
+	  catalog[0].measure[m].dbFlags &= ~ID_MEAS_NOCAL;
 	}
-	
-	catalog[0].measure[m].Mcal = 0;
-	catalog[0].measure[m].dbFlags &= 0xff00;
-	catalog[0].measure[m].dbFlags &= ~ID_MEAS_POOR_PHOTOM;
-	catalog[0].measure[m].dbFlags &= ~ID_MEAS_SKIP_PHOTOM;
-	catalog[0].measure[m].dbFlags &= ~ID_MEAS_AREA;
-	catalog[0].measure[m].dbFlags &= ~ID_MEAS_NOCAL;
       }
     }
@@ -68,5 +77,5 @@
 
   off_t i, k, m;
-  int ecode;
+  int ecode, found, Ns;
   off_t Ntot, Ntry, Nkeep, Nskip;
   float mag;
@@ -95,5 +104,9 @@
       /* skipped via NOCAL, don't mark as skipped */
       ecode = GetPhotcodeEquivCodebyCode (catalog[0].measure[m].photcode);
-      if (ecode != photcode[0].code) continue;
+      found = FALSE;
+      for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
+	if (ecode == photcodes[Ns][0].code) found = TRUE;
+      }
+      if (!found) continue;
 
       /* skip measurements by time range */