Commit 45f8c4e2 authored by LEFEBVREJP email's avatar LEFEBVREJP email
Browse files

Initial working version of f71stream reading.

parent a3a4b543
Pipeline #12654 failed with stages
in 19 minutes and 24 seconds
......@@ -9,5 +9,209 @@ F71Stream::F71Stream(const std::string& file)
: mFile(file)
, mStream(file.c_str())
{
mStream.setReverseBytes(false);
}
class F71Case::PImpl
{
public:
std::vector<int> light_element_nuclides;
std::vector<int> actinide_nuclides;
std::vector<int> fission_product_nuclides;
std::vector<float> light_element_nuclide_abundances;
std::vector<float> actinide_nuclide_abundances;
std::vector<float> fission_product_nuclide_abundances;
std::string title;
std::string basis;
std::vector<float> light_element_gamma_energy_bounds;
std::vector<float> actinide_gamma_energy_bounds;
std::vector<float> fission_product_gamma_energy_bounds;
std::vector<float> light_element_spectrum;
std::vector<float> actinide_spectrum;
std::vector<float> fission_product_spectrum;
std::vector<float> gamma_energy_bounds;
std::vector<float> photon_per_second_spectrum;
std::vector<float> mev_per_second_spectrum;
std::string photon_spectrum_title;
std::vector<float> spontaneous_fission_neutron_source;
std::vector<float> actinide_alpha_n_source;
std::vector<float> total_neutron_spectrum;
std::vector<float> alpha_n_spectrum;
std::vector<float> spontaneous_fission_spectrum;
std::vector<float> neutron_energy_bounds;
double time;
float power;
float flux;
};
F71Case::F71Case()
: p(new PImpl(), [](PImpl* impl) { delete impl; })
{
}
F71Case::F71Case(const F71Case& orig)
: p(new PImpl(), [](PImpl* impl) { delete impl; })
{
p->light_element_nuclides = orig.p->light_element_nuclides;
p->actinide_nuclides = orig.p->actinide_nuclides;
p->fission_product_nuclides = orig.p->fission_product_nuclides;
p->light_element_nuclide_abundances =
orig.p->light_element_nuclide_abundances;
p->actinide_nuclide_abundances = orig.p->actinide_nuclide_abundances;
p->fission_product_nuclide_abundances =
orig.p->fission_product_nuclide_abundances;
p->title = orig.p->title;
p->basis = orig.p->basis;
p->light_element_gamma_energy_bounds =
orig.p->light_element_gamma_energy_bounds;
p->actinide_gamma_energy_bounds = orig.p->actinide_gamma_energy_bounds;
p->fission_product_gamma_energy_bounds =
orig.p->fission_product_gamma_energy_bounds;
p->light_element_spectrum = orig.p->light_element_spectrum;
p->actinide_spectrum = orig.p->actinide_spectrum;
p->fission_product_spectrum = orig.p->fission_product_spectrum;
p->gamma_energy_bounds = orig.p->gamma_energy_bounds;
p->photon_per_second_spectrum = orig.p->photon_per_second_spectrum;
p->mev_per_second_spectrum = orig.p->mev_per_second_spectrum;
p->photon_spectrum_title = orig.p->photon_spectrum_title;
p->spontaneous_fission_neutron_source =
orig.p->spontaneous_fission_neutron_source;
p->actinide_alpha_n_source = orig.p->actinide_alpha_n_source;
p->total_neutron_spectrum = orig.p->total_neutron_spectrum;
p->alpha_n_spectrum = orig.p->alpha_n_spectrum;
p->spontaneous_fission_spectrum = orig.p->spontaneous_fission_spectrum;
p->neutron_energy_bounds = orig.p->neutron_energy_bounds;
p->time = orig.p->time;
p->power = orig.p->power;
p->flux = orig.p->flux;
}
void F71Case::setLightElementNuclides(const std::vector<int>& zaids)
{
p->light_element_nuclides = zaids;
}
void F71Case::setActinideNuclides(const std::vector<int>& zaids)
{
p->actinide_nuclides = zaids;
}
void F71Case::setFissionProductNuclides(const std::vector<int>& zaids)
{
p->fission_product_nuclides = zaids;
}
void F71Case::setLightElementNuclideAbundances(const std::vector<float>& abs)
{
p->light_element_nuclide_abundances = abs;
}
void F71Case::setActinideNuclideAbundances(const std::vector<float>& abs)
{
p->actinide_nuclide_abundances = abs;
}
void F71Case::setFissionProductNuclideAbundances(const std::vector<float>& abs)
{
p->fission_product_nuclide_abundances = abs;
}
void F71Case::setTitle(const std::string& title) { p->title = title; }
void F71Case::setBasis(const std::string& basis) { p->basis = basis; }
void F71Case::setLightElementGammaEnergyBounds(const std::vector<float>& bounds)
{
p->light_element_gamma_energy_bounds = bounds;
}
void F71Case::setActinideGammaEnergyBounds(const std::vector<float>& bounds)
{
p->actinide_gamma_energy_bounds = bounds;
}
void F71Case::setFissionProductGammaEnergyBounds(
const std::vector<float>& bounds)
{
p->fission_product_gamma_energy_bounds = bounds;
}
void F71Case::setLightElementSpectrum(const std::vector<float>& spectrum)
{
p->light_element_spectrum = spectrum;
}
void F71Case::setActinideSpectrum(const std::vector<float>& spectrum)
{
p->actinide_spectrum = spectrum;
}
void F71Case::setFissionProductSpectrum(const std::vector<float>& spectrum)
{
p->fission_product_spectrum = spectrum;
}
void F71Case::setGammaEnergyBounds(const std::vector<float>& bounds)
{
p->gamma_energy_bounds = bounds;
}
void F71Case::setPhotonPerSSpectrum(const std::vector<float>& spectrum)
{
p->photon_per_second_spectrum = spectrum;
}
void F71Case::setMevPerSSpectrum(const std::vector<float>& spectrum)
{
p->mev_per_second_spectrum = spectrum;
}
void F71Case::setPhotonSpectrumTitle(const std::string& title)
{
p->photon_spectrum_title = title;
}
void F71Case::setSpontaneousFissionNeutronSource(
const std::vector<float>& source)
{
p->spontaneous_fission_neutron_source = source;
}
void F71Case::setActinideAlphaNSource(const std::vector<float>& source)
{
p->actinide_alpha_n_source = source;
}
void F71Case::setTotalNeutronSpectrum(const std::vector<float>& spectrum)
{
p->total_neutron_spectrum = spectrum;
}
void F71Case::setAlphaNSpectrum(const std::vector<float>& spectrum)
{
p->alpha_n_spectrum = spectrum;
}
void F71Case::setSpontaneousFissionSpectrum(const std::vector<float>& spectrum)
{
p->spontaneous_fission_spectrum = spectrum;
}
void F71Case::setNeutronEnergyBounds(const std::vector<float>& bounds)
{
p->neutron_energy_bounds = bounds;
}
void F71Case::setTime(double seconds) { p->time = seconds; }
void F71Case::setPower(float power) { p->power = power; }
void F71Case::setFlux(float flux) { p->flux = flux; }
} // namespace radix
......@@ -4,6 +4,7 @@
#include <algorithm>
#include <cmath>
#include <iomanip>
#include <memory>
#include <sstream>
#include <string>
#include <vector>
......@@ -24,10 +25,97 @@ class F71Stream
F71Stream() = delete;
F71Stream(const std::string& mFile);
/**
* @brief read_case Reads a single case from the stream
* @param f71case_type& inventory
*
* Requires the f71case_type to have the following methods:
*
* setLightElementNuclides(const std::vector<int>& zaids)
* setActinideNuclides(const std::vector<int>& zaids)
* setFissionProductNuclides(const std::vector<int>& zaids)
*
* setLightElementNuclideAbundances(const std::vector<float>& abs)
* setActinideNuclideAbundances(const std::vector<float>& abs)
* setFissionProductNuclideAbundances(const std::vector<float>& abs)
*
* setTitle(const std::string& title)
* setBasis(const std::string& basis)
*
* setLightElementGammaEnergyBounds(const std::vector<float>& bounds)
* setActinideGammaEnergyBounds(const std::vector<float>& bounds)
* setFissionProductGammaEnergyBounds(const std::vector<float>& bounds)
*
* setLightElementSpectrum(const std::vector<float>& spectrum)
* setActinideSpectrum(const std::vector<float>& spectrum)
* setFissionProductSpectrum(const std::vector<float>& spectrum)
*
* setGammaEnergyBounds(const std::vector<float>& bounds)
* setPhotonPerSSpectrum(const std::vector<float>& spectrum)
* setMevPerSSpectrum(const std::vector<float>& spectrum)
*
* setPhotonSpectrumTitle(const std::string& title)
* setSpontFissionNeutronSource(const std::vector<float>& source)
* setActinideAlphaNSource(const std::vector<float>& source)
* setTotalNeutronSpectrum(const std::vector<float>& spectrum)
* setAlphaNSpectrum(const std::vector<float>& spectrum)
* setSpontaneousFissionSpectrum(const std::vector<float>& spectrum)
* setNeutronEnergyBounds(const std::vector<float>& bounds)
*
* setTime(double seconds)
* setPower(float power)
* setFlux(float flux)
* @return bool on if inventory was populated
*/
template <typename f71case_type>
bool load_next(f71case_type& inventory);
bool read_case(f71case_type& inventory);
}; // class F71Stream
class F71Case
{
// forward declare private implementation
class PImpl;
// unique pointer to private implmentation
std::unique_ptr<PImpl, void (*)(PImpl*)> p;
public:
F71Case();
F71Case(const F71Case& orig);
void setLightElementNuclides(const std::vector<int>& zaids);
void setActinideNuclides(const std::vector<int>& zaids);
void setFissionProductNuclides(const std::vector<int>& zaids);
void setLightElementNuclideAbundances(const std::vector<float>& abs);
void setActinideNuclideAbundances(const std::vector<float>& abs);
void setFissionProductNuclideAbundances(const std::vector<float>& abs);
void setTitle(const std::string& title);
void setBasis(const std::string& basis);
void setLightElementGammaEnergyBounds(const std::vector<float>& bounds);
void setActinideGammaEnergyBounds(const std::vector<float>& bounds);
void setFissionProductGammaEnergyBounds(const std::vector<float>& bounds);
void setLightElementSpectrum(const std::vector<float>& spectrum);
void setActinideSpectrum(const std::vector<float>& spectrum);
void setFissionProductSpectrum(const std::vector<float>& spectrum);
void setGammaEnergyBounds(const std::vector<float>& bounds);
void setPhotonPerSSpectrum(const std::vector<float>& spectrum);
void setMevPerSSpectrum(const std::vector<float>& spectrum);
void setPhotonSpectrumTitle(const std::string& title);
void setSpontaneousFissionNeutronSource(const std::vector<float>& source);
void setActinideAlphaNSource(const std::vector<float>& source);
void setTotalNeutronSpectrum(const std::vector<float>& spectrum);
void setAlphaNSpectrum(const std::vector<float>& spectrum);
void setSpontaneousFissionSpectrum(const std::vector<float>& spectrum);
void setNeutronEnergyBounds(const std::vector<float>& bounds);
void setTime(double seconds);
void setPower(float power);
void setFlux(float flux);
}; // F71Case
} // namespace radix
#include "radixio/f71stream.i.hh"
......
......@@ -11,8 +11,81 @@
namespace radix
{
template <typename f71case_type>
bool F71Stream::load_next(f71case_type& inventory)
bool F71Stream::read_case(f71case_type& inventory)
{
//
// Load the inventory with processed data
// ITXX - total number of nuclides or ‘0’to signal an end-of-file
// ILITE - number of nuclides in light element library - 0 if not using
// IACT - number of nuclides in actinide library - 0 if not using
// IFP - number of nuclides in fission product library - 0 if not using
// NRFLAG - flag that determines position numbers for data set to be saved
// MSUB1 - time step number of the data to be saved
// NSTEPO - position of requested data in concentration file
// KASEPO - position of the requested case in the concentration file
// JOBPOS - position in the file of the job containing the requested data
// NOCS - subcase number in the case containing the requested data
// NOBLND - 0/1/N – job termination/no blending/N streams blended
// NDSET - nuclear data library unit number
// NTYPE - the type of library on NDSET, 0/>0 card-image / binary library
// NGRP - the number of neutron energy groups
// NELEM - maximum atomic number of nuclides in library
// NVERT - −N/0/N – request integral option and writing N/not wanted/ reading
// N NG -the number of photon energy groups in the master photon library MMN
// - number of irradiation time-step intervals in a subcase MOUT - controls
// writing concentrations on unit NXTR
// If MOUT=0 and MSUB=m, concentration for time-step m is written.
// If MOUT= −N and MSUB>= 0, concentrations for time-steps MSUB
// through N, inclusive, are written.
// INDEX - 0/1 – power read in 58* array/ flux read in 59* array
// MSTAR - time step at which concentration print cutoffs are applied
// NUNIT - 1/2/3/4/5/6 – time (in 60*) units s/min/h/d/y/other
// KBLEND - time-step number indicating when blending is to occur
// NENLE - number of groups in light element photon library
// NENAC - number of groups in actinide photon library
// NENFP - number of groups in fission product photon library
// L1 - max0(1, NENLE)
// L2 - max0(1, NENAC)
// L3 - max0(1, NANFP)
// L4 - max0(1, IACT)
// TMO - time at start of first time-step
// FRACPW - ratio of subcase irradiation time to total case irradiation time
// TCONST - If NUNIT=6, number of seconds in the special time unit.
// TUNIT - time unit printed in output (4-character variable)
// TWRITE - value of time-step accumulative time of concentrations
// written in file
// PWRITE - value of time-step power of concentrations written in file
// FWRITE - value of time-step flux of concentrations written in the file
// The second record type description follows:
// NUCL(ITXX) - ITXX values of nuclide ID’s as ZAS
// (Z=atomic number, A=mass number, S=metastable state, S=0 for ground state)
// X(ITXX) - ITXX values of gram-atom concentrations, at step requested MSUB.
// TITLE(20) - 80 character title for this subcase if NTI=1 or 3
// BASIS(10) - 40 character title for basis of calculation
// - may be blank if NTI=2 or 3
// EGAMLE(NENLE+1)- Gamma group boundaries for light elements
// EGAMAC(NENAC+1)- Gamma group boundaries for actinides
// EGAMFP(NENFP +1)- Gamma group boundaries for fission products
// SPECLE(L1) - L1 values of gamma source spectrum for
// light elements in photons/s
// SPECAC(L2) - L2 values of gamma source spectrum for actinides in photons/s
// SPECFP(L3) - L3 values of gamma source spectrum for
// fission products in photons/s
// ENER(NG+1) - NG +1
// values of gamma energy group boundaries on master photon library
// DSAV(NG) - NG values of master photon spectrum in photons/s
// ESAV(NG) - NG values of master photon spectrum in MeV/s
// ITSAV(20) - 80-character photon spectrum title
// SPNNUC(L4) - L4 values of spontaneous fission neutron source in neutrons/s
// ALPNUC(L4) - L4 values of actinide nuclide (alpha,n) source in neutrons/s
// SPNEUT(NGRP) - NGRP values of total neutron spectrum in neutrons/s
// SPECAN(NGRP) - NGRP values of (alpha,n) spectrum in neutrons/s
// SPECSP(NGRP) - NGRP values of spontaneous fission neutron
// spectrum in neutrons/s
// ENEUTS(NGRP+1) - Neutron energy-group boundaries
// check if the file is open
if (!mStream.is_open())
{
......@@ -27,7 +100,9 @@ bool F71Stream::load_next(f71case_type& inventory)
// header
mStream.readHeader();
// read record content
int num_nuclides = mStream.readInt(); // ITXX
int num_nuclides = mStream.readInt(); // ITXX
if (num_nuclides == 0) return false; // EOF
int num_light_nuclides = mStream.readInt(); // ILITE
int num_actinides = mStream.readInt(); // IACT =
int num_fission_products = mStream.readInt(); // IFP =
......@@ -48,7 +123,7 @@ bool F71Stream::load_next(f71case_type& inventory)
mStream.readInt(); // MOUT =
mStream.readInt(); // INDEX =
mStream.readInt(); // MSTAR =
mStream.readInt(); // NUNIT =
int nunit = mStream.readInt(); // NUNIT =
mStream.readInt(); // KBLEND =
int num_light_element_groups = mStream.readInt(); // NENLE =
int num_actinide_groups = mStream.readInt(); // NENAC =
......@@ -57,27 +132,40 @@ bool F71Stream::load_next(f71case_type& inventory)
int l2 = mStream.readInt(); // L2 =
int l3 = mStream.readInt(); // L3 =
int l4 = mStream.readInt(); // L4 =
mStream.readFloat(); // TMO =
mStream.readFloat(); // FRACPW =
mStream.readFloat(); // TCONST =
mStream.readString(4); // TUNIT =
mStream.readFloat(); // TWRITE =
mStream.readFloat(); // PWRITE =
mStream.readFloat(); // FWRITE =
mStream.readFloat(); // TMO =
mStream.readFloat(); // FRACPW =
mStream.readFloat(); // TCONST =
mStream.readString(4); // TUNIT =
float twrite = mStream.readFloat(); // TWRITE =
float pwrite = mStream.readFloat(); // PWRITE =
float fwrite = mStream.readFloat(); // FWRITE =
// read end of record
mStream.readFooter();
radix_line("Total nuclides " << num_nuclides << "=[" << num_light_nuclides
<< ", " << num_actinides << ", "
<< num_fission_products << "]");
radix_line("Time: " << twrite << " " << nunit << " power(" << pwrite
<< ") flux(" << fwrite << ")");
// read record 2: the nuclides and their amounts
// read header
mStream.readHeader();
std::vector<int> iza_list((size_t(num_nuclides)));
mStream >> iza_list; // values
std::vector<float> abundances((size_t(num_nuclides)));
mStream >> abundances; // gram-atom concentrations
// break iza_list into light/actinide/fission product nuclides lists
std::vector<int> light_nuclides((size_t(num_light_nuclides)));
std::vector<int> actinide_nuclides((size_t(num_actinides)));
std::vector<int> fission_product_nuclides((size_t(num_fission_products)));
mStream >> light_nuclides >> actinide_nuclides >> fission_product_nuclides;
// break abundances by light/actinide/fission product
std::vector<float> light_abs((size_t(num_light_nuclides)));
std::vector<float> actinide_abs((size_t(num_actinides)));
std::vector<float> fission_product_abs((size_t(num_fission_products)));
mStream >> light_abs >> actinide_abs >> fission_product_abs;
std::string title = trim_string(mStream.readString(80));
std::string basis = trim_string(mStream.readString(40));
radix_line("Title(" << title << ") Basis(" << basis << ")");
std::vector<float> gamma_le_bounds((size_t(num_light_element_groups)));
mStream >> gamma_le_bounds;
......@@ -136,7 +224,52 @@ bool F71Stream::load_next(f71case_type& inventory)
mStream.peek(); // peek into the next bit. This ensures EOF is seen in good()
// call
inventory.setLightElementNuclides(light_nuclides);
inventory.setLightElementNuclideAbundances(light_abs);
inventory.setActinideNuclides(actinide_nuclides);
inventory.setActinideNuclideAbundances(actinide_abs);
inventory.setFissionProductNuclides(fission_product_nuclides);
inventory.setFissionProductNuclideAbundances(fission_product_abs);
switch (nunit)
{
case 1: // seconds
inventory.setTime(twrite);
break;
case 2: // minutes
inventory.setTime(twrite * 60.);
break;
case 3: // hours
inventory.setTime(twrite * 3600.);
break;
case 4: // days
inventory.setTime(twrite * 86400);
break;
case 5: // years
inventory.setTime(twrite * 3.154e+7);
break;
}
inventory.setPower(pwrite);
inventory.setFlux(fwrite);
inventory.setTitle(title);
inventory.setBasis(basis);
inventory.setLightElementGammaEnergyBounds(gamma_le_bounds);
inventory.setActinideGammaEnergyBounds(gamma_ae_bounds);
inventory.setFissionProductGammaEnergyBounds(gamma_fpe_bounds);
inventory.setLightElementSpectrum(spectrum_le);
inventory.setActinideSpectrum(spectrum_ac);
inventory.setFissionProductSpectrum(spectrum_fp);
inventory.setGammaEnergyBounds(energy_bounds);
inventory.setPhotonPerSSpectrum(dsav);
inventory.setMevPerSSpectrum(esav);
inventory.setPhotonSpectrumTitle(itsav);
inventory.setSpontaneousFissionNeutronSource(spnnuc);
inventory.setActinideAlphaNSource(alpnuc);
inventory.setTotalNeutronSpectrum(spneut);
inventory.setAlphaNSpectrum(specan);
inventory.setSpontaneousFissionSpectrum(specsp);
inventory.setNeutronEnergyBounds(eneuts);
}
} // namespace radix
#endif /** RADIX_RADIXIO_F71STREAM_I_HH_ */
\ No newline at end of file
#endif /** RADIX_RADIXIO_F71STREAM_I_HH_ */
......@@ -3,4 +3,5 @@ INCLUDE(GoogleTest)
ADD_GOOGLE_TEST(tstCSVFile.cc NP 1)
ADD_GOOGLE_TEST(tstCFGFile.cc NP 1)
ADD_GOOGLE_TEST(tstEafstream.cc NP 1)
ADD_GOOGLE_TEST(tstF71Stream.cc NP 1)
ADD_GOOGLE_TEST(tstHysplitCDump.cc NP 1)
#include "gtest/gtest.h"
#include "radixcore/system.hh"
#include "radixio/f71stream.hh"
#include <fstream>
using namespace radix;
TEST(Radixio, DISABLED_F71Stream)
{
F71Stream stream(
std::string("/home/jap/projects/build/daem/Debug/daemmodel/tests/"
"delfic.Site001.f71"));
std::vector<F71Case> cases;
bool result = true;
int num_cases = 0;
while (result)
{
F71Case single_case;
result = stream.read_case(single_case);
cases.emplace_back(single_case);
num_cases++;
}
std::cout << "Read " << num_cases << " f71 cases." << std::endl;
}
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