From 33f3e8a24d752e7c9494d30493275995d2931ade Mon Sep 17 00:00:00 2001 From: Andreas Adelmann Date: Thu, 11 Aug 2011 15:01:16 +0000 Subject: [PATCH] more tools must later go to OPAL --- .gitattributes | 2 + tools/h5PartDcToVtk.cc | 653 ++++++++++++++++++++++++++++++++++++ tools/h5PartSurfaceToVtk.cc | 532 +++++++++++++++++++++++++++++ 3 files changed, 1187 insertions(+) create mode 100755 tools/h5PartDcToVtk.cc create mode 100644 tools/h5PartSurfaceToVtk.cc diff --git a/.gitattributes b/.gitattributes index 0a776fc..18c6ccb 100644 --- a/.gitattributes +++ b/.gitattributes @@ -560,6 +560,8 @@ tools/H5PartMerge/src/optparse.cpp -text tools/H5PartMerge/src/optparse.hh -text tools/Makefile.am -text tools/README -text +tools/h5PartDcToVtk.cc -text +tools/h5PartSurfaceToVtk.cc -text tools/h5hutcc.in -text tools/h5pAttrib.cc -text tools/h5pToGNUplot.cc -text diff --git a/tools/h5PartDcToVtk.cc b/tools/h5PartDcToVtk.cc new file mode 100755 index 0000000..3329364 --- /dev/null +++ b/tools/h5PartDcToVtk.cc @@ -0,0 +1,653 @@ +/* h5ToVtk.cc + Andreas Adelmann + +*/ + +#include +#include "H5hut.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +//#include +//#include +using namespace std; + +#define MAX_LEN 100 + +/* Function headers */ +int get_option(int argc, const char **argv, const char *opts, const struct long_options *l_opts); +static void print_help(); +static void variable_assign(int argc, const char *argv[]); + +/* Global variables */ +static char* input_name = NULL; +static char* output_name = NULL; +static bool flg_alive = false; +static double z_pos = 0.0; +static int print_all = 0; + +/* `get_option' variables */ +int opt_err = 1; /*get_option prints errors if this is on */ +int opt_ind = 1; /*token pointer */ +const char *opt_arg = NULL; /*flag argument (or value) */ + +/* indication whether the flag (option) requires an argument or not */ +enum { + no_arg = 0, /* doesn't take an argument */ + require_arg, /* requires an argument */ +}; + +/* struct for flags (options) */ +typedef struct long_options +{ + const char *name; /* name of the long option */ + int has_arg; /* whether we should look for an arg */ + char shortval; /* the shortname equivalent of long arg + * this gets returned from get_option */ +} long_options; + +/* List of options in single characters */ +static const char *s_opts = "h1:2:i:o:a:"; + +/* List of options in full words */ +static struct long_options l_opts[] = + { + { "help", no_arg, 'h' }, // Print help page + { "input", require_arg, 'i' }, // Takes input file name + { "output", require_arg, 'o' }, // Takes output file name (without this flag, the program will print to stdout) + { "alive ", no_arg , 'a' }, // also generate the alive dark current to display the dark current source + { NULL, 0, '\0' } + }; + + + +/************************************************************************************ + *********************************** FUNCTIONS ************************************* + *************************************************************************************/ + + +string convert2Int(int number) { + stringstream ss; + ss << setw(5) << setfill('0') << number; + return ss.str(); +} + + +/* get_option is the parsing function that was majorly ported from h5dump utility */ +int get_option(int argc, const char **argv, const char *opts, const struct long_options *l_opts) { + static int sp = 1; /* character index in current token */ + int opt_opt = '?'; /* option character passed back to user */ + + if (sp == 1) + { + /* check for more flag-like tokens */ + if (opt_ind >= argc || argv[opt_ind][0] != '-' || argv[opt_ind][1] == '\0') + { + return EOF; + } + else if (strcmp(argv[opt_ind], "--") == 0) + { + opt_ind++; + return EOF; + } + } + + if (sp == 1 && argv[opt_ind][0] == '-' && argv[opt_ind][1] == '-') + { + /* long command line option */ + const char *arg = &argv[opt_ind][2]; + int i; + + for (i = 0; l_opts && l_opts[i].name; i++) + { + size_t len = strlen(l_opts[i].name); + + if (strncmp(arg, l_opts[i].name, len) == 0) + { + /* we've found a matching long command line flag */ + opt_opt = l_opts[i].shortval; + + if (l_opts[i].has_arg != no_arg) + { + if (arg[len] == '=') + { + opt_arg = &arg[len + 1]; + } + else if (opt_ind < (argc - 1) && argv[opt_ind + 1][0] != '-') + { + opt_arg = argv[++opt_ind]; + } + else if (l_opts[i].has_arg == require_arg) + { + if (opt_err) + fprintf(stderr, "%s: option required for \"--%s\" flag\n", argv[0], arg); + + opt_opt = '?'; + } + } + else + { + if (arg[len] == '=') + { + if (opt_err) + fprintf(stderr, "%s: no option required for \"%s\" flag\n", argv[0], arg); + + opt_opt = '?'; + } + + opt_arg = NULL; + } + + break; + } + } + + if (l_opts[i].name == NULL) + { + /* exhausted all of the l_opts we have and still didn't match */ + if (opt_err) + fprintf(stderr, "%s: unknown option \"%s\"\n", argv[0], arg); + + opt_opt = '?'; + } + + opt_ind++; + sp = 1; + } + else + { + register char *cp; /* pointer into current token */ + + /* short command line option */ + opt_opt = argv[opt_ind][sp]; + + if (opt_opt == ':' || (cp = strchr(opts, opt_opt)) == 0) + { + + if (opt_err) + fprintf(stderr, "%s: unknown option \"%c\"\n", argv[0], opt_opt); + /* if no chars left in this token, move to next token */ + if (argv[opt_ind][++sp] == '\0') + { + opt_ind++; + sp = 1; + } + + return '?'; + } + + if (*++cp == ':') + { + + /* if a value is expected, get it */ + if (argv[opt_ind][sp + 1] != '\0') + { + /* flag value is rest of current token */ + opt_arg = &argv[opt_ind++][sp + 1]; + } + else if (++opt_ind >= argc) + { + if (opt_err) + { + fprintf(stderr, "%s: value expected for option \"%c\"\n", argv[0], opt_opt); + } + opt_opt = '?'; + } + else + { + /* flag value is next token */ + opt_arg = argv[opt_ind++]; + } + + sp = 1; + } + else + { + /* set up to look at next char in token, next time */ + if (argv[opt_ind][++sp] == '\0') + { + /* no more in current token, so setup next token */ + opt_ind++; + sp = 1; + } + + opt_arg = NULL; + } + } + + /* return the current flag character found */ + return opt_opt; +} + +/* Assigns functions according to the parsed result */ +static void variable_assign(int argc, const char *argv[]) +{ + int option; + + /* set options according to the command line */ + while ((option = get_option(argc, argv, s_opts, l_opts)) != EOF) + { + switch ((char)option) + { + case 'h': // Print help page + print_help(); + exit(1); + case 'o': // Print number of steps + output_name = strdup(opt_arg); + break; + case 'i': // Print shorter version without the values + input_name = strdup(opt_arg); + break; + case 'a': // + { flg_alive = true; + z_pos = atof(strdup(opt_arg)); + } + break; + default: + print_help(); + exit(1); + } + } +} + +/* For printing help page */ +static void print_help() +{ + fflush(stdout); + fprintf(stdout, "\nusage: h5ToVtk -i INPUTFILE -o OUTPUTFILE [OPTIONAL_FLAGS] -a ZVALUE\n"); + fprintf(stdout, "\n"); + fprintf(stdout, " FLAGS\n"); + fprintf(stdout, " -h, --help Print help page\n"); + fprintf(stdout, " -i file, --input file (REQUIRED) Takes input base file name to \"file\" (extension h5 is assumed \n"); + fprintf(stdout, " -o file, --output file (REQUIRED) Takes output base file name to \"file\" (extension vtk is added)\n"); + fprintf(stdout, " -a zvalue Only display particles which have servived and reached z value \n"); + + fprintf(stdout, "\n"); + fprintf(stdout, " Examples:\n"); + fprintf(stdout, "\n"); + fprintf(stdout, " /h5ToVtk -i ctf3-injector-darkcurrent-1 -o ctf3-injector-darkcurrent-1- \n"); + fprintf(stdout, "\n"); + fprintf(stdout, " /h5ToVtk -i ctf3-injector-darkcurrent-1 -o ctf3-injector-darkcurrent-1- -a 0.19 \n"); + fprintf(stdout, "\n"); +} + + +int main(int argc, const char *argv[]) +{ + + h5_file_t *h5file = NULL; + + std::ofstream of, ofalive, ofenergy, ofpnum; + + int j; + + int num_dataset; + + int ntime_step = 0; + + variable_assign(argc, argv); + + if(input_name == NULL) { + fprintf(stdout, "missing input file name\n"); + print_help(); + exit(1); + } + + if(output_name == NULL) { + fprintf(stdout, "missing output file name\n"); + print_help(); + exit(1); + } + + string ifn = string(input_name) + string(".h5"); + + h5file = H5OpenFile(ifn.c_str(), H5_O_WRONLY); + + if( h5file == NULL ) { + fprintf(stdout, "unable to open file %s\n", input_name); + print_help(); + exit(1); + } + + ntime_step = H5GetNumSteps(h5file); + + string ffenergy = string(output_name) + string("energy") + string(".dat");//new + ofenergy.open(ffenergy.c_str());//new + assert(ofenergy.is_open());//new + ofenergy << setprecision(10)//new + << "# Energy" << endl;//new + string ffpnum = string(output_name) + string("pnum") + string(".dat");//new + ofpnum.open(ffpnum.c_str());//new + assert(ofpnum.is_open());//new + ofpnum << setprecision(10)//new + << "# time" <<" partNum"<< endl;//new + + if (flg_alive) { + + set idSet; + + + H5SetStep(h5file, ntime_step-1); + num_dataset = H5GetNumDatasets(h5file); + h5_int64_t nparticles = H5GetNumParticles(h5file); + + h5_int64_t* larray = (h5_int64_t*)malloc(sizeof(h5_int64_t)*nparticles); + H5PartReadDataInt64(h5file, "id", larray); + + h5_float64_t* z = (h5_float64_t*)malloc(sizeof(h5_float64_t)*nparticles); + H5PartReadDataFloat64(h5file, "z", z); + + for(unsigned long int n = 0; n < nparticles; ++n) { + if (z[n] >= z_pos) + idSet.insert(larray[n]); + } + + cout << "Last timestep contains " << nparticles << " particles" << endl; + + for (size_t j = 0; j x_alive; + + h5_float64_t* y = (h5_float64_t*)malloc(sizeof(h5_float64_t)*nparticles); + H5PartReadDataFloat64(h5file, "y", y); + vector y_alive; + + h5_float64_t* z = (h5_float64_t*)malloc(sizeof(h5_float64_t)*nparticles); + H5PartReadDataFloat64(h5file, "z", z); + vector z_alive; + + h5_int64_t* ptype = (h5_int64_t*)malloc(sizeof(h5_int64_t)*nparticles); + H5PartReadDataInt64(h5file, "ptype", ptype); + vector ptype_alive; + + for (size_t i = 0; i < nparticles ; i ++) { + if ( idSet.find(larray[i]) != idSet.end() ) { + x_alive.push_back(x[i]); + y_alive.push_back(y[i]); + z_alive.push_back(z[i]); + ptype_alive.push_back(ptype[i]); + } + } + cout<<" ptype_alive size = "<< ptype_alive.size()<Esmall;//smaller than 28eV + vectorEmultip;//Energy zone of SEY>1. + vectorElarge;//larger than 2100eV + double sey_num;//sey total number + int remained_num;// + for (size_t j=0; j idSet; + H5SetStep(h5file,j); + // char s_name[64]="ENERGY";//new + //h5_float64_t Etmp=0.0;//new + //int var_readE=H5PartReadStepAttrib(h5file,s_name,&Etmp);//new + //num_dataset = H5PartGetNumDatasets(h5file);//new + num_dataset = H5GetNumDatasets(h5file); + h5_int64_t nparticles = H5GetNumParticles(h5file); + //Etmp = (Etmp - 0.51099906)*1000000.0;//new eV + //ofenergy<= 2) { + of << 1.0 << endl; + } + else { + cout<< "Step " <=0 "<< endl; + } + } + + of << "LOOKUP_TABLE " << "mytable " << 3 << endl ; + of << 1.0 <<" "<< 0.0 <<" "<< 0.0 <<" "<< 1.0 << endl; + of << 0.0 <<" "<< 1.0 <<" "<< 0.0 <<" "<< 1.0 << endl; + of << 0.0 <<" "<< 0.0 <<" "<< 1.0 <<" "<< 1.0 << endl; + of << endl; + of.close(); + + h5_int64_t* larray = (h5_int64_t*)malloc(sizeof(h5_int64_t)*nparticles); + H5PartReadDataInt64(h5file, "id", larray); + + h5_float64_t* larrayPx = (h5_float64_t*)malloc(sizeof(h5_float64_t)*nparticles); + H5PartReadDataFloat64(h5file, "px", larrayPx); + + h5_float64_t* larrayPy = (h5_float64_t*)malloc(sizeof(h5_float64_t)*nparticles); + H5PartReadDataFloat64(h5file, "py", larrayPy); + + h5_float64_t* larrayPz = (h5_float64_t*)malloc(sizeof(h5_float64_t)*nparticles); + H5PartReadDataFloat64(h5file, "pz", larrayPz); + if(j!=(ntime_step-1)) { + H5PartSetStep(h5file,j+1); + int num_dataset_temp = H5GetNumDatasets(h5file); + h5_int64_t nparticles_temp = H5GetNumParticles(h5file); + + h5_int64_t* larray_temp = (h5_int64_t*)malloc(sizeof(h5_int64_t)*nparticles_temp); + H5PartReadDataInt64(h5file, "id", larray_temp); + + + + for (size_t i = 0; i < nparticles_temp ; i ++) { + idSet.insert(larray_temp[i]); + } + for (size_t i = 0; i < nparticles ; i ++) { + + if(idSet.find(larray[i])==idSet.end()) { + double Etemp=0.51099906*1000000.0*(sqrt(1.0+larrayPx[i]*larrayPx[i]+larrayPy[i]*larrayPy[i]+larrayPz[i]*larrayPz[i])-1); + + double sey_temp=2.469*exp(-0.0004644*Etemp)-1.9*exp(-0.01054*Etemp); + sey_num+=sey_temp; + if ( std::isinf(Etemp)){ + + cout<<"Etemp in time step: "< +#include +#include +#include +#include +#include +#include +#include +//#include +//#include +using namespace std; + +#define MAX_LEN 100 + +/* Function headers */ +int get_option(int argc, const char **argv, const char *opts, const struct long_options *l_opts); +static void print_help(); +static void variable_assign(int argc, const char *argv[]); + +/* Global variables */ +static char* input_name = NULL; +static char* mytemp = NULL; +static char* geo_name = NULL; +static char* output_name = NULL; +static bool flg_alive = false; +static double z_pos = 0.0; +static int print_all = 0; +static int step_ini = 1; +static int dump_freq = 1; +/* `get_option' variables */ +int opt_err = 1; /*get_option prints errors if this is on */ +int opt_ind = 1; /*token pointer */ +const char *opt_arg = NULL; /*flag argument (or value) */ + +/* indication whether the flag (option) requires an argument or not */ +enum { + no_arg = 0, /* doesn't take an argument */ + require_arg, /* requires an argument */ +}; + +/* struct for flags (options) */ +typedef struct long_options +{ + const char *name; /* name of the long option */ + int has_arg; /* whether we should look for an arg */ + char shortval; /* the shortname equivalent of long arg + * this gets returned from get_option */ +} long_options; + +/* List of options in single characters */ +static const char *s_opts = "h1:2:i:g:o:f:d:"; + +/* List of options in full words */ +static struct long_options l_opts[] = + { + { "help", no_arg, 'h' }, // Print help page + { "input", require_arg, 'i' }, // Takes input file name + { "geometry", require_arg, 'g' }, // Takes input geometry file name + { "output", require_arg, 'o' }, // Takes output file name (without this flag, the program will print to stdout) + { "first", require_arg, 'f' }, // first step + { "dumpfreq", require_arg, 'd' }, // dump frequency of the + { NULL, 0, '\0' } + }; + + + +/************************************************************************************ + *********************************** FUNCTIONS ************************************* + *************************************************************************************/ + + +string convert2Int(int number) { + stringstream ss; + ss << setw(5) << setfill('0') << number; + return ss.str(); +} + + +/* get_option is the parsing function that was majorly ported from h5dump utility */ +int get_option(int argc, const char **argv, const char *opts, const struct long_options *l_opts) { + static int sp = 1; /* character index in current token */ + int opt_opt = '?'; /* option character passed back to user */ + + if (sp == 1) + { + /* check for more flag-like tokens */ + if (opt_ind >= argc || argv[opt_ind][0] != '-' || argv[opt_ind][1] == '\0') + { + return EOF; + } + else if (strcmp(argv[opt_ind], "--") == 0) + { + opt_ind++; + return EOF; + } + } + + if (sp == 1 && argv[opt_ind][0] == '-' && argv[opt_ind][1] == '-') + { + /* long command line option */ + const char *arg = &argv[opt_ind][2]; + int i; + + for (i = 0; l_opts && l_opts[i].name; i++) + { + size_t len = strlen(l_opts[i].name); + + if (strncmp(arg, l_opts[i].name, len) == 0) + { + /* we've found a matching long command line flag */ + opt_opt = l_opts[i].shortval; + + if (l_opts[i].has_arg != no_arg) + { + if (arg[len] == '=') + { + opt_arg = &arg[len + 1]; + } + else if (opt_ind < (argc - 1) && argv[opt_ind + 1][0] != '-') + { + opt_arg = argv[++opt_ind]; + } + else if (l_opts[i].has_arg == require_arg) + { + if (opt_err) + fprintf(stderr, "%s: option required for \"--%s\" flag\n", argv[0], arg); + + opt_opt = '?'; + } + } + else + { + if (arg[len] == '=') + { + if (opt_err) + fprintf(stderr, "%s: no option required for \"%s\" flag\n", argv[0], arg); + + opt_opt = '?'; + } + + opt_arg = NULL; + } + + break; + } + } + + if (l_opts[i].name == NULL) + { + /* exhausted all of the l_opts we have and still didn't match */ + if (opt_err) + fprintf(stderr, "%s: unknown option \"%s\"\n", argv[0], arg); + + opt_opt = '?'; + } + + opt_ind++; + sp = 1; + } + else + { + register char *cp; /* pointer into current token */ + + /* short command line option */ + opt_opt = argv[opt_ind][sp]; + + if (opt_opt == ':' || (cp = strchr(opts, opt_opt)) == 0) + { + + if (opt_err) + fprintf(stderr, "%s: unknown option \"%c\"\n", argv[0], opt_opt); + /* if no chars left in this token, move to next token */ + if (argv[opt_ind][++sp] == '\0') + { + opt_ind++; + sp = 1; + } + + return '?'; + } + + if (*++cp == ':') + { + + /* if a value is expected, get it */ + if (argv[opt_ind][sp + 1] != '\0') + { + /* flag value is rest of current token */ + opt_arg = &argv[opt_ind++][sp + 1]; + } + else if (++opt_ind >= argc) + { + if (opt_err) + { + fprintf(stderr, "%s: value expected for option \"%c\"\n", argv[0], opt_opt); + } + opt_opt = '?'; + } + else + { + /* flag value is next token */ + opt_arg = argv[opt_ind++]; + } + + sp = 1; + } + else + { + /* set up to look at next char in token, next time */ + if (argv[opt_ind][++sp] == '\0') + { + /* no more in current token, so setup next token */ + opt_ind++; + sp = 1; + } + + opt_arg = NULL; + } + } + + /* return the current flag character found */ + return opt_opt; +} + +/* Assigns functions according to the parsed result */ +static void variable_assign(int argc, const char *argv[]) +{ + int option; + + /* set options according to the command line */ + while ((option = get_option(argc, argv, s_opts, l_opts)) != EOF) + { + switch ((char)option) + { + case 'h': // Print help page + print_help(); + exit(1); + case 'o': // Print number of steps + output_name = strdup(opt_arg); + break; + case 'i': // Print shorter version without the values + input_name = strdup(opt_arg); + break; + case 'g': // + { geo_name = strdup(opt_arg); + + } + break; + case 'f': // + { + + //fflush(stdout); + //char ss* = NULL; + //ss = strdup(opt_arg); + mytemp = strdup(opt_arg); + step_ini = (int)strtod(mytemp,NULL); + fprintf(stdout,"string=%s, %s \n",opt_arg, mytemp); + fprintf(stdout,"step_ini=%i \n",step_ini); + + } + break; + case 'd': // + { + + //fflush(stdout); + //char df* = NULL; + //df = strdup(opt_arg); + fprintf(stdout,"string=%s \n",opt_arg); + mytemp = strdup(opt_arg); + dump_freq = (int)strtod(mytemp,NULL); + fprintf(stdout,"string=%s \n",opt_arg); + fprintf(stdout,"dump_freq=%i \n",dump_freq); + + } + break; + default: + print_help(); + exit(1); + } + } +} + +/* For printing help page */ +static void print_help() +{ + fflush(stdout); + fprintf(stdout, "\nusage: h5SurfaceVtk -i INPUTFILE -g GEOMETRYFILE -o OUTPUTFILE\n"); + fprintf(stdout, "\n"); + fprintf(stdout, " FLAGS\n"); + fprintf(stdout, " -h, --help Print help page\n"); + fprintf(stdout, " -i file, --input file (REQUIRED) Takes input base file name to \"file\" (extension h5 is assumed \n"); + fprintf(stdout, " -g file, --input geometry file (REQUIRED) Takes input base file name to \"file\" (extension h5 is assumed \n"); + fprintf(stdout, " -o file, --output file (REQUIRED) Takes output base file name to \"file\" (extension vtk is added)\n"); + + fprintf(stdout, "\n"); + fprintf(stdout, " Examples first dump step 100 and dump frequency 100:\n"); + fprintf(stdout, "\n"); + fprintf(stdout, " /h5SurfaceVtk -i Cavity_Mag_Surface -g ../10 -o p- -f 100 -d 100 \n"); + fprintf(stdout, "\n"); + +} + + +int main(int argc, const char *argv[]) +{ + H5PartFile *h5file = NULL; + + std::ofstream of, ofalive, ofenergy, ofpnum; + + int j; + + int num_dataset; + + int ntime_step = 0; + + variable_assign(argc, argv); + + if(input_name == NULL) { + fprintf(stdout, "missing input file name\n"); + print_help(); + exit(1); + } + + if(output_name == NULL) { + fprintf(stdout, "missing output file name\n"); + print_help(); + exit(1); + } + + string ifn = string(input_name) + string(".h5"); + + h5file = H5PartOpenFile(ifn.c_str(), H5PART_READ); + + if( h5file == NULL ) { + fprintf(stdout, "unable to open file %s\n", input_name); + print_help(); + exit(1); + } + + string gfn = string(geo_name) + string(".h5"); + hid_t plist_id = H5Pcreate(H5P_FILE_ACCESS); //Property list identifier + // H5Pset_fapl_mpio(plist_id, MPI_COMM_WORLD, MPI_INFO_NULL); + // Create a new file collectively and release property list identifier. + hid_t file_id = H5Fopen(gfn.c_str(), H5F_ACC_RDONLY, plist_id); + assert(file_id >= 0); + H5Pclose(plist_id); + int numbfaces_global_m; + ///////////////////////////////////////////// + // Read dataset "surface" from .h5 file + //////////////////////////////////////////// + + hsize_t dimsf[2];//dataset dimensions + herr_t status; + + hid_t dset_id = H5Dopen(file_id, "/surface"); + assert(dset_id >= 0); + // Create the dataspace for the dataset. + hid_t x = H5Dget_space(dset_id); + H5Sget_simple_extent_dims(x, dimsf, NULL); + hid_t filespace = H5Screate_simple(2, dimsf, NULL); + + numbfaces_global_m = dimsf[0]; + int *allbfaces_m = (int*)malloc(numbfaces_global_m * 4*sizeof(int)); + + int *bfaces_idx_m = (int*)malloc(numbfaces_global_m*sizeof(int)); + + //Tribarycent_m = new Vector_t[numbfaces_global_m]; + // Create property list for collective dataset write. + hid_t plist_id2 = H5Pcreate(H5P_DATASET_XFER); + //H5Pset_dxpl_mpio(plist_id2, H5FD_MPIO_COLLECTIVE); + status = H5Dread(dset_id, H5T_NATIVE_INT, H5S_ALL, filespace, plist_id2, allbfaces_m); + assert(status >= 0); + // Store local boundary faces, discard others + int nof_sym_m = 0; // Count number of symmetry planes. + int const nogeosym_flag = 0; // heronion outputs a index, case 0 indicates no symmetry plane, 1 for (x,y) sym, 2 for (y,z), 3 for (z,x),none 0 means the current triangle is not on a real surface of geometry but on a symmetry plane. + + + for(int i = 0; i < numbfaces_global_m; i ++) { + bfaces_idx_m[i] = i; + if(allbfaces_m[4 * i] > nogeosym_flag) { + nof_sym_m += 1; + if(i < numbfaces_global_m - 1) { + for(int j = 0; j < numbfaces_global_m - i; j ++) { + allbfaces_m[4*(i+j)] = allbfaces_m[4*(i+j+1)]; + allbfaces_m[4*(i+j)+1] = allbfaces_m[4*(i+j+1)+1]; + allbfaces_m[4*(i+j)+2] = allbfaces_m[4*(i+j+1)+2]; + allbfaces_m[4*(i+j)+3] = allbfaces_m[4*(i+j+1)+3]; + } + } else + numbfaces_global_m = numbfaces_global_m - 1; + } + } + H5Dclose(dset_id); + H5Sclose(filespace); + //////////////////////////////// + // Also read dataset "coords" + //////////////////////////////// + hsize_t dimsf_c[2]; + herr_t status_c; + hid_t dset_id_c = H5Dopen(file_id, "/coords"); + assert(dset_id_c >= 0); + + // Create the dataspace for the dataset. + hid_t cp_space = H5Dget_space(dset_id_c); + H5Sget_simple_extent_dims(cp_space, dimsf_c, NULL); + hid_t filespace_c = H5Screate_simple(2, dimsf_c, NULL); + + int numpoints_global_m = dimsf_c[0]; + double *point_coords = (double*)malloc(3 * numpoints_global_m*sizeof(double)); + hid_t plist_id3 = H5Pcreate(H5P_DATASET_XFER); + //H5Pset_dxpl_mpio(plist_id3, H5FD_MPIO_COLLECTIVE); + status_c = H5Dread(dset_id_c, H5T_NATIVE_DOUBLE, H5S_ALL, filespace_c, plist_id3, point_coords); + assert(status_c >= 0); + + + + + double **geo3Dcoords_m;//= malloc(sizeof(double)*numpoints_global_m*3); + geo3Dcoords_m = (double**)malloc(sizeof(double)*numpoints_global_m); + for (int i=0;i idSet; + size_t step_local = step_ini + (j-1)*dump_freq; + H5PartSetStep(h5file,step_local); + num_dataset = H5PartGetNumDatasets(h5file); + h5part_int64_t triNum = H5PartGetNumParticles(h5file); + double qi = 0.0; + H5PartReadStepAttrib ( h5file,"qi",&qi); + cout << "Working on timestep " << step_local << endl; + + h5part_float64_t* PrimaryLoss = (h5part_float64_t*)malloc(sizeof(h5part_float64_t)*triNum); + H5PartReadDataFloat64(h5file, "PrimaryLoss", PrimaryLoss); + + h5part_float64_t* SecondaryLoss = (h5part_float64_t*)malloc(sizeof(h5part_float64_t)*triNum); + H5PartReadDataFloat64(h5file, "SecondaryLoss", SecondaryLoss); + + h5part_float64_t* FNEmissionLoss = (h5part_float64_t*)malloc(sizeof(h5part_float64_t)*triNum); + H5PartReadDataFloat64(h5file, "FNEmissionLoss", FNEmissionLoss); + + h5part_int64_t* TriID = (h5part_int64_t*)malloc(sizeof(h5part_int64_t)*triNum); + H5PartReadDataInt64(h5file, "TriangleID", TriID); + + string ffn = string("vtk/") + string(output_name) + string("Surface-") + convert2Int(j) + string(".vtk"); + + of.open(ffn.c_str()); + assert(of.is_open()); + of.precision(6); + of << "# vtk DataFile Version 2.0" << endl; // first line of a vtk file + of << "generated using DataSink::writeGeoToVtk" << endl; // header + of << "ASCII" << endl << endl; // file format + of << "DATASET UNSTRUCTURED_GRID" << endl; // dataset structrue + of << "POINTS " << numpoints_global_m << " float" << endl; // data num and type + for(int i = 0; i < numpoints_global_m ; i ++) + of << geo3Dcoords_m[i][0] << " " << geo3Dcoords_m[i][1] << " " << geo3Dcoords_m[i][2] << endl; + of << endl; + + of << "CELLS " << numbfaces_global_m << " " << 4 * numbfaces_global_m << endl; + for(int i = 0; i < numbfaces_global_m ; i ++) + of << "3 " << allbfaces_m[4*i+1] << " " << allbfaces_m[4*i+2] << " " << allbfaces_m[4*i+3] << endl; + of << "CELL_TYPES " << numbfaces_global_m << endl; + for(int i = 0; i < numbfaces_global_m ; i ++) + of << "5" << endl; + of << "CELL_DATA " << numbfaces_global_m << endl; + of << "SCALARS " << "PrimaryLoss" << " float " << "1" << endl; + of << "LOOKUP_TABLE " << "default" << endl ; + for(int i = 0; i < numbfaces_global_m ; i ++) + of << fabs(PrimaryLoss[i]/qi) << endl; + of << "SCALARS " << "SecondaryLoss" << " float " << "1" << endl; + of << "LOOKUP_TABLE " << "default" << endl ; + for(int i = 0; i < numbfaces_global_m ; i ++) + of << fabs(SecondaryLoss[i]/qi) << endl; + of << "SCALARS " << "FNEmissionLoss" << " float " << "1" << endl; + of << "LOOKUP_TABLE " << "default" << endl ; + for(int i = 0; i < numbfaces_global_m ; i ++) + of << fabs(FNEmissionLoss[i]/qi) << endl; + of << "SCALARS " << "TotalLoss" << " float " << "1" << endl; + of << "LOOKUP_TABLE " << "default" << endl ; + for(int i = 0; i < numbfaces_global_m ; i ++) + of << fabs((FNEmissionLoss[i]+SecondaryLoss[i]+PrimaryLoss[i])/qi) << endl; + of << endl; + + of.close(); + + //cout <<"ptype size = "<< sizeof(ptype) << endl; + + free(PrimaryLoss); + free(SecondaryLoss); + free(FNEmissionLoss); + free(TriID); + } + free(point_coords); + for (int i=0;i