diff --git a/sammy/src/io/CMakeLists.txt b/sammy/src/io/CMakeLists.txt index 43ab2c968fb2177fe4ae41c02e6004c9305a8346..d5af49a65b26d7d5553f8d94f530ca1a67561397 100644 --- a/sammy/src/io/CMakeLists.txt +++ b/sammy/src/io/CMakeLists.txt @@ -2,20 +2,25 @@ INCLUDE(SammyPackageSetup) TRIBITS_PACKAGE(io) +find_package(HDF5 REQUIRED COMPONENTS C CXX HL) configure_file(${CMAKE_CURRENT_SOURCE_DIR}/../sammy_conf.h.in sammy_conf.h @ONLY) -INCLUDE_DIRECTORIES(${CMAKE_CURRENT_BINARY_DIR}) - +INCLUDE_DIRECTORIES(${CMAKE_CURRENT_BINARY_DIR} + ${HDF5_INCLUDE_DIRS}) # from find_package() # Set headers SET(HEADERS ExperimentalParametersIO.h + OdfIO.h + Hdf5IO.h interface/cpp/ExperimentalParametersIOInterface.h ) APPEND_SET(IO_SOURCES ExperimentalParametersIO.cpp + OdfIO.cpp + Hdf5IO.cpp interface/cpp/ExperimentalParametersIOInterface.cpp interface/fortran/ExperimentalParametersIO_M.f90 interface/fortran/ExperimentalParametersIO_I.f90 @@ -34,3 +39,15 @@ INSTALL(FILES ${HEADERS} DESTINATION "src/io") TRIBITS_ADD_TEST_DIRECTORIES(tests) TRIBITS_PACKAGE_POSTPROCESS() + +target_link_libraries(SammyIOUtilsLib ${HDF5_LIBRARIES} ) # from find_package() + +# -------------------------------------------- +# Build separate executable to convert ODF +# to HDF5 files, install to normal location +# -------------------------------------------- + +add_executable(convertodf convertodf.cpp) +target_link_libraries(convertodf SammyIOUtilsLib) +install(PROGRAMS ${CMAKE_CURRENT_BINARY_DIR}/convertodf DESTINATION "bin/") + diff --git a/sammy/src/io/Hdf5IO.cpp b/sammy/src/io/Hdf5IO.cpp new file mode 100644 index 0000000000000000000000000000000000000000..b81cc466aee0927a331f377c9230ffe978e4bb1e --- /dev/null +++ b/sammy/src/io/Hdf5IO.cpp @@ -0,0 +1,104 @@ +#include "Hdf5IO.h" + +namespace sammy { + + Hdf5IO::Hdf5IO() {} + + Hdf5IO::Hdf5IO(std::map<std::string,int> &header){ + // copy function var directly to class header member + fillHDF5header(header); + } + + void Hdf5IO::writeODFtoHDF5(std::map<std::string,int> odfhead, + std::vector<std::vector<double>> &data, + std::string hdf5Filename){ + + // ---------------------------------------------------------- + // Set up data sections + // ---------------------------------------------------------- + if( data.size() != 11 ){ + std::cerr << "---------------------------------------" << std::endl; + std::cerr << "Reading of ODF files with < 11 sections" + << "has not been implemented yet. " << std::endl; + std::cerr << "---------------------------------------" << std::endl; + throw std::runtime_error("Unrecognized ODF file data."); + } + + // attempt to make std::strings + + // HDF5 requires either char * or char[] + std::vector<const char *> sect11 = {"energy","trans","dtrans", + "cro","dcro","cs","zero","co","zero","idk","dcs"}; + std::vector<const char *> sect9; + + // HDF5 accepts set-size double arrays for datasets: copy over + const size_t ns = odfhead["numSection"]; + const size_t nc = odfhead["numChan"]; + double dset_data[ns][nc]; + for( size_t i=0;i<data.size();++i ){ + std::copy(data[i].begin(),data[i].end(), dset_data[i]); + } + + // ---------------------------------------------------------- + // prepare 3 datasets: data, data section names, and header + // ---------------------------------------------------------- + + // --- data space ----- + hsize_t ddim[2]; + ddim[0] = odfhead["numSection"]; + ddim[1] = odfhead["numChan"]; + int drank = sizeof(ddim) / sizeof(hsize_t); + std::string dataname("/data/sections"); + + // --- section names space ----- + hsize_t sdim[1]; + sdim[0] = odfhead["numSection"]; + int srank = sizeof(sdim) / sizeof(hsize_t); + std::string sectionnames("/data/sectnames"); + + // --- header ----- + fillHDF5header(odfhead); + + hsize_t mdim[1]; + mdim[0] = h5head.size(); + int mrank = sizeof(mdim) / sizeof(hsize_t); + + const std::string headername("/meta/header"); + const std::string headKeyword("key"); + const std::string headVal("value"); + + // --- define compound type for header struct ----- + H5::CompType htype(sizeof(sammy::headKeypair)); + htype.insertMember(headKeyword, HOFFSET(sammy::headKeypair, key ),H5::StrType(H5::PredType::C_S1, H5T_VARIABLE)); + htype.insertMember(headVal, HOFFSET(sammy::headKeypair, value ),H5::PredType::NATIVE_INT); + + // ---------------------------------------------------------- + // open file and write + // ---------------------------------------------------------- + H5::DataSpace dspace(drank,ddim); + H5::DataSpace sspace(srank,sdim); + H5::DataSpace mspace(mrank,mdim); + H5::H5File file = H5::H5File(hdf5Filename,H5F_ACC_TRUNC); + H5::Group group = H5::Group( file.createGroup( "/data" )); + H5::Group metaGroup = H5::Group( file.createGroup( "/meta" )); + + H5::DataSet dataset = H5::DataSet(file.createDataSet(dataname,H5::PredType::NATIVE_DOUBLE,dspace)); + H5::DataSet nameset = H5::DataSet(file.createDataSet(sectionnames,H5::StrType(H5::PredType::C_S1, H5T_VARIABLE),sspace)); + H5::DataSet metaset = H5::DataSet(file.createDataSet(headername,htype,mspace)); + + dataset.write( dset_data, H5::PredType::NATIVE_DOUBLE ); + nameset.write( sect11.data(), H5::StrType(H5::PredType::C_S1, H5T_VARIABLE) ); + metaset.write( h5head.data(), htype ); + } + + void Hdf5IO::fillHDF5header(std::map<std::string,int> &odfhead){ + + for( auto & member : odfhead ){ + headKeypair hd; + hd.key = member.first.c_str(); + hd.value = odfhead[member.first]; + h5head.push_back(hd); + } + } + +} \ No newline at end of file diff --git a/sammy/src/io/Hdf5IO.h b/sammy/src/io/Hdf5IO.h new file mode 100644 index 0000000000000000000000000000000000000000..dd073a0f4a1b8785c94721d30ff260a961390bd0 --- /dev/null +++ b/sammy/src/io/Hdf5IO.h @@ -0,0 +1,60 @@ +#ifndef SAMMY_HDF5IO_H +#define SAMMY_HDF5IO_H + +#include <vector> +#include <fstream> +#include <iostream> +#include <cstring> +#include <stdexcept> +#include "OdfIO.h" +#include "H5Cpp.h" + +namespace sammy{ + + /***************************************************************** + * This struct defines the keypair values that will be used to + * compose the header for HDF5 metadata + ****************************************************************/ + struct headKeypair{ + const char * key; + int value; + }; + + /***************************************************************** + * The Hdf5IO class is designed to house any HDF5 input/output + * routines inlcuding functions to convert ORELA Data Format (ODF) + * files into HDF5 files. + ****************************************************************/ + class Hdf5IO { + + public: + Hdf5IO(); + Hdf5IO(std::map<std::string,int> &header); + virtual ~Hdf5IO() {} + + /************************************************************* + * Write data from ODF file into HDF5 + * + * @param head the header information from the ODF file + * @param data the sections of data read from the ODF file + * @param hdf5Filename name of HDF5 file to be written + ************************************************************/ + void writeODFtoHDF5(std::map<std::string,int> odfhead, + std::vector<std::vector<double>> &data, + std::string hdf5Filename); + + /************************************************************* + * Move data into an HDF5 friendly data-structure: array of + * keyword and value structs + * + * @param head the header information from the ODF file + * @param h5head the array to be filled and written to HDF5 + ************************************************************/ + void fillHDF5header(std::map<std::string,int> &odfhead); + + private: + std::vector<headKeypair> h5head; + }; +} + +#endif /* SAMMY_HDF5IO_H */ \ No newline at end of file diff --git a/sammy/src/io/OdfIO.cpp b/sammy/src/io/OdfIO.cpp new file mode 100644 index 0000000000000000000000000000000000000000..11a2241c5094ce20223b3aa4642d2d047e38b7b0 --- /dev/null +++ b/sammy/src/io/OdfIO.cpp @@ -0,0 +1,180 @@ +#include "OdfIO.h" + +namespace sammy{ + + OdfIO::OdfIO() {} + + double OdfIO::readFloat(std::istream &fs, bool unix){ + + unsigned char tmp1,tmp2,tmp3,tmp4; + int i=0; + float value; + + if(!unix){ // wacky VMS + fs.read( (char*)&tmp1, 1 ); // 1 byte each + fs.read( (char*)&tmp2, 1 ); + fs.read( (char*)&tmp3, 1 ); + fs.read( (char*)&tmp4, 1 ); + if((tmp2 & 0x7f) == 0){ + if(tmp1+tmp3+tmp4 != 0) { + tmp1=tmp2=tmp3=0; + } + tmp2=0; + } + else tmp2--; + + i=(tmp2 << 24) + (tmp1 << 16) + (tmp4 << 8) + tmp3; + } + else{ + fs.read( (char*)&tmp1, 1 ); + fs.read( (char*)&tmp2, 1 ); + fs.read( (char*)&tmp3, 1 ); + fs.read( (char*)&tmp4, 1 ); + i=(tmp4 << 24) + (tmp3 << 16) + (tmp2 << 8 ) +tmp1; + } + // grab int pointer, cast to float, dereference + value = *reinterpret_cast<float *>( &i ); + return (double)value; + } + + int OdfIO::readInteger(std::istream &fs){ + size_t intsize = sizeof(int); + int val; + fs.read( (char*)&val, intsize ); + return val; + } + + std::vector<double> OdfIO::getSection(int start, int readSize, int section, + std::istream &in){ + in.seekg(0,std::ios::beg); + std::vector<double> dataSec; + int block,skipped,channel,max; + size_t intsize = sizeof(int); + + dataSec.resize(readSize); + + block = start/125; + channel = start - block * 125; + block = head["ifb"] + block * head["numSection"]; + block--; + skipped=( block + (section-1) )*512; // skip 512 bytes per block + skipped += 4; // account for control word + skipped += channel*4; // and add channels not to read + in.ignore(skipped); // jump skipped-many bytes, alt: in.seekg(skipped,std::ios::cur); + + max = std::min(readSize,125-channel); + + int read = 0; + for(int i=0;i<max;i++){ + if( head["mode"] == 3 ) dataSec[read] = readFloat(in,head["unix"]); + else if( head["mode"] == 1 ) dataSec[read] = readInteger(in); + read++; + } + in.ignore(intsize*2); // skip over two ints in the block which contain no data + + while( read<readSize ){ + skipped=(head["numSection"]-1)*512; // skip over other sections + in.ignore(skipped); + in.ignore(intsize); // discard control word + + for( int i=0; i<125 && read<readSize; i++ ) { + if( head["mode"] == 3) dataSec[read] = readFloat(in,head["unix"]); + else if(head["mode"] == 1) dataSec[read] = readInteger(in); + read++; + } + in.ignore(intsize*2); // skip over two ints in the block which contain no data + } + + return dataSec; + } + + std::vector<std::vector<double>> OdfIO::getAllSections(std::istream &in){ + // --- go to start of file --- + in.seekg(0, std::ios::beg); + + std::vector<std::vector<double>> data; + data.resize(head["numSection"]); + + // --- grab each section --- + for( int section=1; section<=head["numSection"]; ++section ){ + data[section-1].resize(head["numChan"]); + data[section-1] = getSection(0,head["numChan"],section,in); + } + + return data; + } + + void OdfIO::readHeader(std::istream &in, int unix){ + size_t intsize = sizeof(int); + int bytes_per_block = 512; + int bytes_per_word = 4; + head["unix"] = unix; + + // comments based on OPRODF manual: Jack Craven, report no. ORNL/CSD/TM-45, 1978 + // ** keep in mind that changes may have happened between that report and + // ** the writing of the files!!! For example it seems that 36 bit ints are no + // ** longer used (32 bit now). + in.read( (char*)&head["mode"], intsize ); // 0,1,3 for 18-bit ints, 36-bit ints, floating-point vals + in.read( (char*)&head["source"], intsize ); // 0,1,2 for SEL, CSISRS, or ENDF/B data + in.read( (char*)&head["irun"], intsize ); // run number + in.read( (char*)&head["ncblks"], intsize ); // starting block no. for the comment section + in.read( (char*)&head["ncwrds"], intsize ); // number of words in comment section + in.read( (char*)&head["nsblks"], intsize ); // starting block for scalar section + in.read( (char*)&head["nswrds"], intsize ); // number of words in scalar section + in.read( (char*)&head["ncstrt"], intsize ); // starting word number of DAQ scalar and counter section in scalar section + in.read( (char*)&head["ncntrs"], intsize ); // number of words in DAQ scalar and counter section + in.read( (char*)&head["nxstrt"], intsize ); // starting word number of DAQ variable section in scalar section + in.read( (char*)&head["nxwrds"], intsize ); // number of words in DAQ variable section + in.read( (char*)&head["npblks"], intsize ); // starting block number of parameter section + in.read( (char*)&head["npwrds"], intsize ); // number of words in parameter section + in.read( (char*)&head["ndtype"], intsize ); // 0,1 for data section described by table in parameter section, data section is dependent on dataset 1 (i.e. datasets 2-x correspond to energy in dataset 1) + in.read( (char*)&head["numSection"], intsize ); // number of sections in each dataset + in.read( (char*)&head["ifb"], intsize ); // + in.read( (char*)&head["numChan"], intsize ); // number of channels in each section + in.read( (char*)&head["zan"], intsize ); // ENDF/B charge/mass id + in.read( (char*)&head["awr"], intsize ); // ENDF/B ratio of nuclear mass (iso,ele,molec,mix) to neutron + in.read( (char*)&head["mat"], intsize ); // ENDF/B MAT number + in.read( (char*)&head["mf"], intsize ); // ENDF/B file number + in.read( (char*)&head["mt"], intsize ); // ENDF/B reaction type + in.read( (char*)&head["varswt"], intsize ); // 0,1 for dataset 1 decreasing, and increasing in value (ignored if ndtype is 0) + in.read( (char*)&head["dctswt"], intsize ); // 0,1 for not dead-time corrected, and dead-time corrected + in.read( (char*)&head["strt"], intsize ); // starting word number of each dataset in data section (mode 0 only) + in.read( (char*)&head["ndwend"], intsize ); // ending word number of last word written in block NDBEND + + + /* TODO: the following doesn't work properly yet, we ignore for now */ + // ------------------------------------ + // --- Read crunch table + // ------------------------------------ + // energy index table. each word contains the largest energy for each n blocks in the data set. N=(NDWRDS/125)+1 + ndtbl.resize(100); + for( int i=0;i<100;++i ) in.read( (char*)&ndtbl[i], intsize ); + // ------------------------------------ + // --- comments section --- + // ------------------------------------ + int bytes_to_comments = (head["ncblks"]) * bytes_per_block; + std::vector<unsigned char> comments; + int bytes_in_comments = head["ncwrds"]*bytes_per_word; + comments.resize( bytes_in_comments ); + + in.seekg(bytes_to_comments,std::ios::beg); + in.ignore(bytes_per_word); + + for( int i=0;i<bytes_in_comments;++i ){ + in.read( (char*)&comments[i], 1 ); + comments[i] = comments[i]; + } + // ------------------------------------ + // --- scalars section --- + // ------------------------------------ + int bytes_to_scalars = (head["nsblks"]) * bytes_per_block; + } + + void OdfIO::printHeader(){ + for( auto member : head ){ + std::cout << member.first << " = " << head[member.first] << std::endl; + } + + } + +} \ No newline at end of file diff --git a/sammy/src/io/OdfIO.h b/sammy/src/io/OdfIO.h new file mode 100644 index 0000000000000000000000000000000000000000..9dc497a2e31708a64b730ef8b490829ecec584dc --- /dev/null +++ b/sammy/src/io/OdfIO.h @@ -0,0 +1,89 @@ +#ifndef SAMMY_ODFIO_H +#define SAMMY_ODFIO_H + +#include <vector> +#include <fstream> // ifstream +#include <iostream> +#include <map> + + +namespace sammy{ + + /***************************************************************** + * The OdfIO class houses ORELA Data Format (ODF) reading routines + * so that old data files can be read and used for SAMMY. + ****************************************************************/ + class OdfIO { + + public: + OdfIO(); + virtual ~OdfIO() {} + + /************************************************************* + * Read a single floating point value, based on whether it is in + * VAX-VMS format (?) or unix, from an ODF file + * + * @param fs a reference to the ODF file stream + * @param unix whether or not the format is unix-style + * @return val the stored value + ************************************************************/ + double readFloat(std::istream &fs, bool unix); + + /************************************************************* + * Read a single integer + * @param fs a reference to the ODF file stream + * @return val the stored value + ************************************************************/ + int readInteger(std::istream &fs); + + /************************************************************* + * Read one section of data. Sections are typcially energy, transmission, + * cross section, etc. + * + * @param start what channel to start in + * @param readSize how many channels to read + * @param section which section to read from + * @param head a struct containing the header info + * @param in the input data file stream + * @return dataSec all requested channels of requested section + ************************************************************/ + std::vector<double> getSection(int start, int readSize, int section, + std::istream &in); + + /************************************************************* + * Read all sections of data into double-vector, grabbing all info. + * This calls getSection() for all sections and channels + * + * @param head a struct containing the header info + * @param in the input data file stream + * @return data all channels and all sections of data + ************************************************************/ + std::vector<std::vector<double>> getAllSections(std::istream &in); + + /************************************************************* + * Read the ODF header + * + * @param fs a reference to the ODF file stream + * @param unix 1 or 0 for true or false unix format + ************************************************************/ + void readHeader(std::istream &in, int unix); + + const std::map<std::string, int> getHeader() const{ + return head; + } + + /************************************************************* + * Print the ODF header + * + * @param head an ODF-header struct + ************************************************************/ + void printHeader(); + + private: + std::map<std::string, int> head; + std::vector<int> ndtbl; + + }; +} + +#endif /* SAMMY_ODFIO_H */ \ No newline at end of file diff --git a/sammy/src/io/convertodf.cpp b/sammy/src/io/convertodf.cpp new file mode 100644 index 0000000000000000000000000000000000000000..ea652389e43ae387960079f7713bfc941af5dfa1 --- /dev/null +++ b/sammy/src/io/convertodf.cpp @@ -0,0 +1,54 @@ +#include <iostream> +#include "H5Cpp.h" +#include "OdfIO.h" +#include "Hdf5IO.h" + +int main(int argc, char const *argv[]) +{ + // ------------------------------------------------------------------------- + // --- Grabbing the input args + // ------------------------------------------------------------------------- + if( argc < 3 ){ + std::cerr << "Wrong number of arguments given." << std::endl + << "Specify input ODF file and output HDF5 file." << std::endl; + return -1; + } + + std::string odfFilename = argv[1]; + std::string hdf5Filename = argv[2]; + + std::cout << "Reading from file: " << odfFilename << std::endl; + std::cout << "Writing to file: " << hdf5Filename << std::endl; + + // ------------------------------------------------------------------------- + // --- Read ODF file + // ------------------------------------------------------------------------- + sammy::OdfIO odf; + std::vector<std::vector<double>> data; + + // --- Open file --- + std::ifstream in( odfFilename, std::ios::in | std::ios::binary ); + // --- kill bad open --- + if( !in.is_open() ){ + std::cerr << "Error opening file: " << odfFilename << std::endl; + return 1; + } + + int unix = false; + // --- read header --- + odf.readHeader(in,unix); + + // --- read data sections --- + data = odf.getAllSections(in); + + in.close(); + + // ------------------------------------------------------------------------- + // --- Output HDF5 file + // ------------------------------------------------------------------------- + + sammy::Hdf5IO h5writer; + h5writer.writeODFtoHDF5(odf.getHeader(),data,hdf5Filename); + + return 0; +} \ No newline at end of file diff --git a/sammy/src/io/tests/CMakeLists.txt b/sammy/src/io/tests/CMakeLists.txt index 792c9ba7979c88cae0921c30411b4bceada90d12..2318d6e1225a364c123ef105208160163a043f15 100644 --- a/sammy/src/io/tests/CMakeLists.txt +++ b/sammy/src/io/tests/CMakeLists.txt @@ -1,3 +1,4 @@ include(NemesisTest) -ADD_NEMESIS_TEST(ExperimentalParametersIOTest.cpp TIMEOUT 60000 NP 1 ) +ADD_NEMESIS_TEST(ExperimentalParametersIOTest.cpp TIMEOUT 60000 NP 1 ) +ADD_NEMESIS_TEST(OdfIOTest.cpp TIMEOUT 60000 NP 1 ) diff --git a/sammy/src/io/tests/OdfIOTest.cpp b/sammy/src/io/tests/OdfIOTest.cpp new file mode 100644 index 0000000000000000000000000000000000000000..dccc8c14c65ad4da730788f8920f19c2cde066df --- /dev/null +++ b/sammy/src/io/tests/OdfIOTest.cpp @@ -0,0 +1,82 @@ +#include "Nemesis/gtest/Gtest_Functions.hh" +#include "Nemesis/gtest/nemesis_gtest.hh" +#include "../OdfIO.h" + +/******************************************************************************* + * Test whether we read floats with VAX format properly + ******************************************************************************/ +TEST(OdfIO,readFloatVax) +{ + sammy::OdfIO odf; + bool is_unix = false; + double readVal; + float * dummy; + double gold; + char byte1 = 0b00000001; + char byte2 = 0b00000010; + char byte3 = 0b00000011; + char byte4 = 0b00000100; + // | | | | + int vax = 0b00000001000000010000010000000011; + // byte order 2--, 1, 4, 3 + + std::stringstream in; + in << byte1 << byte2 << byte3 << byte4; + + readVal = odf.readFloat(in,is_unix); + + dummy = (float *) &vax; + gold = *dummy; + + ASSERT_EQ( readVal, gold ); +} + +/******************************************************************************* + * Test whether we read floats with UNIX format properly + ******************************************************************************/ +TEST(OdfIO,readFloatUnix) +{ + sammy::OdfIO odf; + bool is_unix = true; + double readVal; + float * dummy; + double gold; + char byte1 = 0b00000001; + char byte2 = 0b00000010; + char byte3 = 0b00000011; + char byte4 = 0b00000100; + // | | | | + int unix = 0b00000100000000110000001000000001; + + std::stringstream in; + in << byte1 << byte2 << byte3 << byte4; + + readVal = odf.readFloat(in,is_unix); + + dummy = (float *) &unix; + gold = *dummy; + + ASSERT_EQ( readVal, gold ); +} + +/******************************************************************************* + * Test whether we read integers right + ******************************************************************************/ +TEST(OdfIO,readInteger) +{ + sammy::OdfIO odf; + int gold = 16909060; + char byte1 = 0b00000001; + char byte2 = 0b00000010; + char byte3 = 0b00000011; + char byte4 = 0b00000100; + int readVal; + + std::stringstream in; + in << byte4 << byte3 << byte2 << byte1; + + readVal = odf.readInteger(in); + + ASSERT_EQ( readVal, gold ); +} +