Commit 8a7ba99f authored by jclamber's avatar jclamber
Browse files

updates on tipsy2gadget

git-svn-id: http://svn.oamp.fr/repos/uns_projects/trunk@86 f264a43e-d52d-4b82-913b-c2bd5215a18a
parent 7babd297
/*
* Example of fortran call:
*
* integer ndark
* double precision time
* real xx(NPARTMAX),yy(NPARTMAX),zz(NPARTMAX),mass(NPARTMAX)
* real vx(NPARTMAX),vy(NPARTMAX),vz(NPARTMAX)
*
* call readDMtip(ndark,xx(1),yy(1),zz(1),vx(1),vy(1),vz(1),mass(1),time)
*/
// ============================================================================
// Copyright Jean-Charles LAMBERT - 2010-2014
// Centre de donneeS Astrophysiques de Marseille (CeSAM)
// e-mail: Jean-Charles.Lambert@lam.fr
// address: Aix Marseille Universite, CNRS, LAM
// Laboratoire d'Astrophysique de Marseille
// Pole de l'Etoile, site de Chateau-Gombert
// 38, rue Frederic Joliot-Curie
// 13388 Marseille cedex 13 France
// CNRS UMR 7326
// ============================================================================
#define _FILE_OFFSET_BITS 64
#include <stdio.h>
......@@ -21,42 +22,29 @@
#include <nemo.h>
#include <uns.h>
// NEMO parameters
// ------------------------------------------------------------
// Nemo variable
const char * defv[] = { // use `::'string because of 'using namespace std'
"in=???\n tipsy input file ",
"out=???\n gadget2 output file)",
"scale=100\n scaling factor (pos,hsml)",
"scale=1.0\n scaling factor (pos,hsml)",
"verbose=f\n verbose mode ",
"VERSION=1.0\n compiled on <"__DATE__"> JCL ",
NULL
};
const char * usage="Simple converter tipsy to gadget2";
void readtipsystd_(const char * name,std::string outname,float scale)
/*
char *name,int *nbods, int *ngas,int *nstar, int* ndark, int *ndim,double *time,
float *mass, float *xx, float *yy, float *zz,
float *vx, float *vy, float *vz,float *soft,float *pot,float *metals,
float *tform, float* rho,float* temp)
*/
void readtipsystd_(const char * filename,std::string outname,float scale)
{
XDR xdr;
FILE *fp;
int i,j;
char filename[80];
// char filename[80];
struct star_particle sp;
struct dark_particle dp;
struct gas_particle gp;
struct dump head;
for (i=0; i<80; i++){
if (name[i]==32) filename[i]=0;
else
filename[i]=name[i];
}
uns::CunsOut * unsout = new uns::CunsOut(outname,"gadget2",false);
fp = fopen(filename,"r");
......@@ -64,12 +52,12 @@ char *name,int *nbods, int *ngas,int *nstar, int* ndark, int *ndim,double *time,
/* read header */
xdr_header(&xdr,&head);
std::cerr << "Scale = "<<scale <<"\n";
/*printf(" Particles in TIPSY file: %d gas, %d dark, %d star\n",
head.nsph,head.ndark,head.nstar);*/
assert(head.ndim==3);
float time = head.time;
int ndim=head.ndim;
//int ndim=head.ndim;
int nbods = head.nbodies;
int ndark=head.ndark;
int ngas = head.nsph;
......@@ -77,44 +65,49 @@ char *name,int *nbods, int *ngas,int *nstar, int* ndark, int *ndim,double *time,
/*printf(" Time: %g\n",*time);*/
fprintf(stderr,"Time=%f nbods=%d nhalo=%d ngas=%d nstars=%d\n",
time,nbods,ndark,ngas,nstar);
time,nbods,ndark,ngas,nstar);
/* reads data for gas particles */
float * mass=NULL, * xx=NULL, * vv=NULL, * rho=NULL, *hsml=NULL;
bool ok;
// set time
ok=unsout->setData("time",time);
if (ok) ; // remove compiler warning
/* reads data for gas particles */
i = 0;
while (i<head.nsph) {
xdr_gas(&xdr,&gp);
if (!i) {
mass = new float[head.nsph];
xx = new float[head.nsph*3];
vv = new float[head.nsph*3];
rho = new float[head.nsph];
hsml = new float[head.nsph];
}
mass[i]=gp.mass;
xx[i*3+0]=gp.pos[0]*scale;
xx[i*3+1]=gp.pos[1]*scale;
xx[i*3+2]=gp.pos[2]*scale;
vv[i*3+0]=gp.vel[0];
vv[i*3+1]=gp.vel[1];
vv[i*3+2]=gp.vel[2];
rho[i]=gp.rho;
//temp[i]=gp.temp;
hsml[i]=gp.hsmooth*scale;
//metals[i]=gp.metals;
//pot[i]=gp.phi;
xdr_gas(&xdr,&gp);
if (!i) {
mass = new float[head.nsph];
xx = new float[head.nsph*3];
vv = new float[head.nsph*3];
rho = new float[head.nsph];
hsml = new float[head.nsph];
}
mass[i] =gp.mass;
xx[i*3+0]=gp.pos[0]*scale;
xx[i*3+1]=gp.pos[1]*scale;
xx[i*3+2]=gp.pos[2]*scale;
vv[i*3+0]=gp.vel[0];
vv[i*3+1]=gp.vel[1];
vv[i*3+2]=gp.vel[2];
rho[i] =gp.rho;
//temp[i]=gp.temp;
hsml[i]=gp.hsmooth*scale;
//metals[i]=gp.metals;
//pot[i]=gp.phi;
++i;
}
if (head.nsph) {
ok=unsout->setData("gas","pos",head.nsph,xx,false);
ok=unsout->setData("gas","vel",head.nsph,vv,false);
ok=unsout->setData("gas","pos" ,head.nsph,xx ,false);
ok=unsout->setData("gas","vel" ,head.nsph,vv ,false);
ok=unsout->setData("gas","mass",head.nsph,mass,false);
ok=unsout->setData("gas","rho",head.nsph,rho,false);
ok=unsout->setData("gas","rho" ,head.nsph,rho ,false);
ok=unsout->setData("gas","hsml",head.nsph,hsml,false);
delete [] mass;
delete [] xx;
......@@ -124,12 +117,10 @@ char *name,int *nbods, int *ngas,int *nstar, int* ndark, int *ndim,double *time,
}
/* reads data for dark particles */
#if 1
j=0;
i = 0;
while (i<head.ndark) {
xdr_dark(&xdr,&dp);
if (!i) {
mass = new float[head.ndark];
......@@ -149,26 +140,24 @@ char *name,int *nbods, int *ngas,int *nstar, int* ndark, int *ndim,double *time,
++j;
}
if (head.ndark) {
ok=unsout->setData("halo","pos",head.ndark,xx,false);
ok=unsout->setData("halo","vel",head.ndark,vv,false);
if (head.ndark) {
ok=unsout->setData("halo","pos" ,head.ndark,xx ,false);
ok=unsout->setData("halo","vel" ,head.ndark,vv ,false);
ok=unsout->setData("halo","mass",head.ndark,mass,false);
delete [] mass;
delete [] xx;
delete [] vv;
}
}
/* reads data for star particles */
j = 0;
i=0;
while (i<head.nstar) {
xdr_star(&xdr,&sp);
xdr_star(&xdr,&sp);
if (!i) {
mass = new float[head.ndark];
xx = new float[head.ndark*3];
vv = new float[head.ndark*3];
mass = new float[head.nstar];
xx = new float[head.nstar*3];
vv = new float[head.nstar*3];
}
mass[i]=sp.mass;
......@@ -180,19 +169,18 @@ if (head.ndark) {
vv[i*3+2]=sp.vel[2];
++i;
++j;
}
if (head.nstar) {
ok=unsout->setData("stars","pos",head.nstar,xx,false);
ok=unsout->setData("stars","vel",head.nstar,vv,false);
if (head.nstar) {
ok=unsout->setData("stars","pos" ,head.nstar,xx ,false);
ok=unsout->setData("stars","vel" ,head.nstar,vv ,false);
ok=unsout->setData("stars","mass",head.nstar,mass,false);
delete [] mass;
delete [] xx;
delete [] vv;
}
}
#endif
fclose(fp);
//printf("Data for %d particles read from %s \n",i,filename);
......@@ -201,6 +189,8 @@ if (head.nstar) {
}
// ------------------------------------------------------------
// main program
int main(int argc, char ** argv) {
// start NEMO
initparam(const_cast<char**>(argv),const_cast<char**>(defv));
......@@ -209,15 +199,14 @@ int main(int argc, char ** argv) {
// Get input parameters
std::string simname = (getparam ((char *) "in" ));
std::string outname = (getparam ((char *) "out" ));
int scale = (getiparam((char *) "scale" ));
float scale = (getdparam((char *) "scale" ));
std::cerr << "Scale = "<<scale <<"\n";
readtipsystd_(simname.c_str(),outname,scale);
// finish NEMO
// finish NEMO
finiparam();
}
// ------------------------------------------------------------
......
......@@ -31,6 +31,7 @@ const char * defv[] = {
"deltav=0.0,0.0,0.0\n velocity of in1 w.r.t. in2",
"shift=1\n Shift over 1st or 2nd one?",
"zerocm=f\n Centering On Mass after stacking?",
"time=0.0\n set output snapshot time",
"verbose=f\n verbose on/off",
"VERSION=1.O\n compiled on <"__DATE__"> JCL ",
NULL,
......@@ -90,11 +91,18 @@ void addComponent(std::string comp, CunsIn * uns1,CunsIn * uns2,
addArray(comp,"u" ,1,uns1,uns2,unsout,verbose);
addArray(comp,"hsml" ,1,uns1,uns2,unsout,verbose);
addArray(comp,"rho" ,1,uns1,uns2,unsout,verbose);
addArray(comp,"temp" ,1,uns1,uns2,unsout,verbose);
addArray(comp,"metal" ,1,uns1,uns2,unsout,verbose);
}
if (comp == "stars") {
addArray(comp,"age" ,1,uns1,uns2,unsout,verbose);
addArray(comp,"metal" ,1,uns1,uns2,unsout,verbose);
}
}
// ------------------------------------------------------------
// process
void process(CunsIn * uns1,CunsIn * uns2, char * out, char * dr, char * dv, bool com, int shift, bool verbose)
void process(CunsIn * uns1,CunsIn * uns2, char * out, char * dr, char * dv, bool com, int shift, bool verbose, float otime)
{
// convert string to vector
std::vector<float> deltar=CSnaptools::stringToVector<float>(dr,3,0.0);
......@@ -126,6 +134,7 @@ void process(CunsIn * uns1,CunsIn * uns2, char * out, char * dr, char * dv, bool
}
if (com)
std::vector<double> com=unsout->snapshot->moveToCom(); // shift to COM
unsout->snapshot->setData("time",otime); // snapshot time
unsout->snapshot->save(); // save file
}
......@@ -144,7 +153,8 @@ int main(int argc, char ** argv )
bool com = getbparam((char *) "zerocm" );
int shift = getiparam((char *) "shift" );
bool verbose = getbparam((char *) "verbose" );
float otime = getdparam((char *) "time" );
// instantiate a new UNS input object (for reading)
uns::CunsIn * uns1 = new uns::CunsIn(in1,"all","all",verbose);
......@@ -170,7 +180,7 @@ int main(int argc, char ** argv )
std::cerr << "UNS files types are not identical, aborting....\n";
std::exit(1);
}
process(uns1,uns2,out,dr,dv,com,shift,verbose);
process(uns1,uns2,out,dr,dv,com,shift,verbose,otime);
} else {
std::cerr << "Can't read nextFrame ....\n";
}
......
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