Commit 1ac67cc9 authored by jclamber's avatar jclamber
Browse files

more polish on fits

git-svn-id: http://svn.oamp.fr/repos/uns_projects/trunk@121 f264a43e-d52d-4b82-913b-c2bd5215a18a
parent 3e1e9d5e
......@@ -15,8 +15,8 @@
#include <assert.h>
#define _vectmath_h // put this statement to avoid conflict with C++ vector class
namespace nos { // use namespace to avoid conflict between
// typedef char * string from stdinc.h
// and using std::string from ccfits.h
// typedef char * string from stdinc.h
// and using std::string from ccfits.h
typedef char * string;
#include "getparam.h"
}
......@@ -29,26 +29,27 @@ typedef char * string;
// The library is enclosed in a namespace.
using namespace CCfits;
// ------------------------------------------------------------
// -------------`-----------------------------------------------
// Nemo variable
const char * defv[] = { // use `::'string because of 'using namespace std'
"in=???\n tipsy input file ",
"out=???\n output file)",
"in=???\n Fits input file ",
"out=???\n output file",
"type=nemo\n nemo|gadget2",
"dmin=\n minimum data value",
"dmax=\n maximum data value",
"zmin=-1\n #z plane min",
"zmax=-1\n #z plane max",
"zmin=\n #z plane min",
"zmax=\n #z plane max",
"verbose=f\n verbose mode",
"VERSION=1.0\n compiled on <"__DATE__"> JCL ",
NULL
};
const char * usage="Simple converter fits to uns compatible output format";
int readImage(std::string in, std::string out, std::string type,int zmin, int zmax);
int readImage(std::string in, std::string out, std::string type);
float dmin = std::numeric_limits<float>::min();
float dmax = std::numeric_limits<float>::max();
int zmin = std::numeric_limits<int>::min();
int zmax = std::numeric_limits<int>::max();
int main(int argc, char ** argv )
{
......@@ -58,8 +59,12 @@ int main(int argc, char ** argv )
// Get input parameters
std::string simname (nos::getparam ((char *) "in" ));
std::string outname (nos::getparam ((char *) "out" ));
int zmin (nos::getiparam ((char *) "zmin" ));
int zmax (nos::getiparam ((char *) "zmax" ));
if (nos::hasvalue((char *) "zmin")) {
zmin=(nos::getiparam ((char *) "zmin" ));
}
if (nos::hasvalue((char *) "zmax")) {
zmax=(nos::getiparam ((char *) "zmax" ));
}
if (nos::hasvalue((char *) "dmin")) {
dmin=(nos::getdparam ((char *) "dmin" ));
}
......@@ -74,24 +79,21 @@ int main(int argc, char ** argv )
else
FITS::setVerboseMode(false);
try {
if (!readImage(simname,outname,typeout,zmin,zmax)) std::cerr << " readImage() \n";
try {
if (!readImage(simname,outname,typeout))
std::cerr << " readImage() ok \n";
}
catch (FitsException&)
// will catch all exceptions thrown by CCfits, including errors
// found by cfitsio (status != 0)
{
catch (FitsException&) {
// will catch all exceptions thrown by CCfits, including errors
// found by cfitsio (status != 0)
std::cerr << " Fits Exception Thrown by test function \n";
}
return 0;
nos::finiparam();
}
int readImage(std::string in, std::string out, std::string type, int zmin, int zmax)
int readImage(std::string in, std::string out, std::string type)
{
// std::auto_ptr<FITS> pInfile(new FITS("/windows7/JDL/ngc6503.fits",Read,true));
std::auto_ptr<FITS> pInfile(new FITS(in,Read,false));
......@@ -102,27 +104,21 @@ int readImage(std::string in, std::string out, std::string type, int zmin, int z
std::cerr << " Begin \n";
std::valarray<float> contents;
#if 1
try {
image = &pInfile->extension(0);
std::cerr << "There is an extension\n";
std::cerr << "#axis ="<<image->axes() << std::endl;
//(pInfile->extension(1)).read(contents);
((ExtHDU *) image)->read(contents);
image = &pInfile->extension(1);
std::cerr << "There is an extension\n";
std::cerr << "#axis ="<<image->axes() << std::endl;
//(pInfile->extension(1)).read(contents);
((ExtHDU *) image)->read(contents);
} catch (FitsException&) {
#endif
image = &pInfile->pHDU();
std::cerr << "There is a pHDU\n";
std::cerr << "#axis ="<<image->axes() << std::endl;
//(pInfile->pHDU()).read(contents);
((PHDU *) image)->read(contents);
image = &pInfile->pHDU();
std::cerr << "There is a pHDU\n";
std::cerr << "#axis ="<<image->axes() << std::endl;
//(pInfile->pHDU()).read(contents);
((PHDU *) image)->read(contents);
}
#if 1
// read all user-specifed, coordinate, and checksum keys in the image
image->readAllKeys();
......@@ -146,31 +142,29 @@ int readImage(std::string in, std::string out, std::string type, int zmin, int z
} else {
az=0;
}
if (zmin==-1) zmin=0;
if (zmax==-1) zmax=az;
std::vector <float> pos,hsml,rho;
long nan=0;
std::cerr << "value = " << contents.size() << "\n";
for (unsigned long i=0; i<contents.size(); i++) {
if (std::isfinite(contents[i]) && contents[i]>dmin && contents[i] < dmax) {
if (std::isfinite(contents[i]) && contents[i]>=dmin && contents[i]<=dmax) {
long z_i=int(i/(ax*ay)); // current Z plane
assert(z_i<=az);
if (z_i>=zmin && z_i<=zmax) { // inside Z selection
long nxy=i-z_i*(ax*ay);// #pixels (x/y) of the latest Z plane
long y_i=int(nxy/ax); // current Y coordinate
assert(y_i<ay);
long x_i=nxy-(y_i*ax); // current X coordinate
long z_i=int(i/(ax*ay)); // current Z plane
assert(z_i<=az);
pos.push_back(x_i*1.0);
pos.push_back(y_i*1.0);
pos.push_back(z_i*1.0);
rho.push_back(contents[i]);
hsml.push_back(0.8);
}
//std::cerr << x_i << " " << y_i << " " << z_i << " " << contents[i] << "\n";
if (z_i>=zmin && z_i<=zmax) { // inside Z selection
long nxy=i-z_i*(ax*ay);// #pixels (x/y) of the latest Z plane
long y_i=int(nxy/ax); // current Y coordinate
assert(y_i<ay);
long x_i=nxy-(y_i*ax); // current X coordinate
assert(z_i<=az);
pos.push_back(x_i*1.0);
pos.push_back(y_i*1.0);
pos.push_back(z_i*1.0);
rho.push_back(contents[i]);
hsml.push_back(0.8);
}
//std::cerr << x_i << " " << y_i << " " << z_i << " " << contents[i] << "\n";
} else { // NAN or INF
nan++;
}
......@@ -187,7 +181,6 @@ int readImage(std::string in, std::string out, std::string type, int zmin, int z
// save snapshot
unsout->snapshot->save();
return 0;
#endif
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment