new tools for handling SDSS mask

mytools/contour_pixels.cpp

@@ -0,0 +1,47 @@
+#include <vector>
+#include <healpix_map.h>
+#include <healpix_map_fitsio.h>
+#include "contour_pixels.hpp"
+
+using namespace std;
+
+static const bool DEBUG = true;
+
+void computeContourPixels(Healpix_Map<float>& m, vector<int> contour)
+{
+  for (int p = 0; p < m.Npix(); p++)
+    {
+      fix_arr<int, 8> result;
+
+      m.neighbors(p, result);
+      for (int q = 0; q < 8; q++)
+	{
+	  if (result[q] < 0)
+	    continue;
+
+	  delta = (m[p]-0.5)*(m[result[q]]-0.5);
+	  if (delta < 0)
+	    {
+	      contour.push_back(p);
+	      // This is boundary go to next pixel
+	      break;
+	    }
+	}      
+    }
+
+  if (DEBUG)
+    {
+      Healpix_Map<int> contour_map;
+
+      contour_map.SetNside(RING, m.Nside());
+      contour_map.fill(0);
+      for (int p = 0; p < contour.size(); p++)
+	{
+	  contour_map[contour[p]]=1;
+	}
+
+      fitshandle h;
+      h.create("!contour_map.fits");
+      write_Healpix_map_to_fits(h, contour_map, planckType<int>());
+    }
+}

mytools/contour_pixels.hpp

@@ -0,0 +1,9 @@
+#ifndef  __CONTOUR_PIXELS_HPP
+#define  __CONTOUR_PIXELS_HPP
+
+#include <vector>
+#include <healpix_map.h>
+
+void computeContourPixels(Healpix_Map<float>& map, std::vector<int> contour);
+
+#endif

mytools/generateFromCatalog.cpp

@@ -0,0 +1,79 @@
+#include <healpix_map.h>
+#include <healpix_map_fitsio.h>
+#include <boost/format.hpp>
+#include <iostream>
+#include <fstream>
+#include <vector>
+#include <string>
+#include "generateFromCatalog_conf.h"
+#include "contour_pixels.hpp"
+
+using namespace std;
+using boost::format;
+
+struct NYU_Data
+{
+  int index;
+  int sector;
+  int region;
+  double ra, dec;
+  double cz;
+  double fgotten;
+  double phi_z;
+};
+
+void loadData(const string& fname, vector<NYU_Data> & data)
+{
+  ifstream f(fname.c_str());
+
+  while (!f.eof())
+    {
+      NYU_Data d;
+      f >> d.index >> d.sector >> d.region >> d.ra >> d.dec >> d.cz >> d.fgotten >> d.phi_z;
+      data.push_back(d);
+    }
+}
+
+
+int main(int argc, char **argv)
+{
+  generateFromCatalog_info args_info;
+  generateFromCatalog_conf_params args_params;
+ 
+  generateFromCatalog_conf_init(&args_info);
+  generateFromCatalog_conf_params_init(&args_params);
+  
+  args_params.check_required = 0;
+  if (generateFromCatalog_conf_ext (argc, argv, &args_info, &args_params))
+    return 1;
+  
+  if (!args_info.configFile_given)
+    {
+      if (generateFromCatalog_conf_required (&args_info, GENERATEFROMCATALOG_CONF_PACKAGE))
+        return 1;
+    }
+  else
+    {
+      args_params.check_required = 1;
+      args_params.initialize = 0;
+      if (generateFromCatalog_conf_config_file (args_info.configFile_arg,
+					 &args_info,
+					 &args_params))
+	return 1;
+    }
+  
+  generateFromCatalog_conf_print_version();
+
+  cout << "Loading NYU data..." << endl;
+  vector<NYU_Data> data;
+  Healpix_Map<float> mask;
+  vector<int> pixel_list;
+  loadData(args_info.catalog_arg, data);
+
+  cout << "Loading mask..." << endl;
+  read_Healpix_map_from_fits(args_info.mask_arg, mask);
+
+  computeContourPixels(mask,pixel_list);
+
+  return 0;
+}

mytools/generateFromCatalog.ggo

@@ -0,0 +1,11 @@
+package "generateFromCatalog"
+version "alpha"
+
+option "configFile" - "Configuration filename" string optional
+
+option "catalog" - "Input NYU-VAGC catalog" string required
+option "mask" - "Healpix mask of unobserved data (in Equatorial coordinates)" string required
+option "density_fake" - "Number density of boundary fake tracers (1 h^3/ Mpc^3)"  double optional default="1"
+
+option "output"	- "Filename of particle datafile" string required
+option "params" - "Output parameters of the datacube" string required