diff --git a/mytools/contour_pixels.cpp b/mytools/contour_pixels.cpp new file mode 100644 index 0000000..725b78c --- /dev/null +++ b/mytools/contour_pixels.cpp @@ -0,0 +1,47 @@ +#include +#include +#include +#include "contour_pixels.hpp" + +using namespace std; + +static const bool DEBUG = true; + +void computeContourPixels(Healpix_Map& m, vector contour) +{ + for (int p = 0; p < m.Npix(); p++) + { + fix_arr 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 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()); + } +} diff --git a/mytools/contour_pixels.hpp b/mytools/contour_pixels.hpp new file mode 100644 index 0000000..341e409 --- /dev/null +++ b/mytools/contour_pixels.hpp @@ -0,0 +1,9 @@ +#ifndef __CONTOUR_PIXELS_HPP +#define __CONTOUR_PIXELS_HPP + +#include +#include + +void computeContourPixels(Healpix_Map& map, std::vector contour); + +#endif diff --git a/mytools/generateFromCatalog.cpp b/mytools/generateFromCatalog.cpp new file mode 100644 index 0000000..7b41380 --- /dev/null +++ b/mytools/generateFromCatalog.cpp @@ -0,0 +1,79 @@ +#include +#include +#include +#include +#include +#include +#include +#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 & 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 data; + Healpix_Map mask; + vector 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; +} diff --git a/mytools/generateFromCatalog.ggo b/mytools/generateFromCatalog.ggo new file mode 100644 index 0000000..bda567d --- /dev/null +++ b/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 \ No newline at end of file