IAP GITLAB

Skip to content
Snippets Groups Projects
Commit 6d650059 authored by Remy Prechelt's avatar Remy Prechelt
Browse files

Restructure output for single run per library.

parent 95688419
No related branches found
No related tags found
1 merge request!317Output infrastructure and Python analysis library.
......@@ -33,7 +33,7 @@ namespace corsika {
setNodes(); // put each particle on stack in correct environment volume
// start this event (i.e. this shower)
output_.startOfEvent();
output_.startOfShower();
while (!stack_.isEmpty()) {
while (!stack_.isEmpty()) {
......@@ -57,7 +57,7 @@ namespace corsika {
}
// end this event (i.e. this shower)
output_.endOfEvent();
output_.endOfShower();
}
template <typename TTracking, typename TProcessList, typename TOutput, typename TStack,
......
......@@ -13,7 +13,9 @@
namespace corsika {
template <typename TOutput>
ObservationPlane<TOutput>::ObservationPlane(Plane const& obsPlane, DirectionVector const& x_axis, bool deleteOnHit)
ObservationPlane<TOutput>::ObservationPlane(Plane const& obsPlane,
DirectionVector const& x_axis,
bool deleteOnHit)
: plane_(obsPlane)
, deleteOnHit_(deleteOnHit)
......@@ -41,9 +43,7 @@ namespace corsika {
auto const displacement = trajectory.getPosition(1) - plane_.getCenter();
// add our particles to the output file stream
this->write(particle.getPID(),
energy,
displacement.dot(xAxis_),
this->write(particle.getPID(), energy, displacement.dot(xAxis_),
displacement.dot(yAxis_),
(trajectory.getPosition(1) - plane_.getCenter()).getNorm());
......@@ -166,46 +166,46 @@ namespace corsika {
template <typename TOutput>
YAML::Node ObservationPlane<TOutput>::getConfig() const {
using namespace units::si;
using namespace units::si;
// construct the top-level node
YAML::Node node;
// construct the top-level node
YAML::Node node;
// basic info
node["type"] = "ObservationPlane";
// basic info
node["type"] = "ObservationPlane";
// the center of the plane
auto const center{plane_.getCenter()};
// the center of the plane
auto const center{plane_.getCenter()};
// save each component in its native coordinate system
auto const center_coords{center.getCoordinates(center.getCoordinateSystem())};
node["plane"]["center"].push_back(center_coords.getX() / 1_m);
node["plane"]["center"].push_back(center_coords.getY() / 1_m);
node["plane"]["center"].push_back(center_coords.getZ() / 1_m);
node["plane"]["center.units"] = "m";
// save each component in its native coordinate system
auto const center_coords{center.getCoordinates(center.getCoordinateSystem())};
node["plane"]["center"].push_back(center_coords.getX() / 1_m);
node["plane"]["center"].push_back(center_coords.getY() / 1_m);
node["plane"]["center"].push_back(center_coords.getZ() / 1_m);
node["plane"]["center.units"] = "m";
// the normal vector of the plane
auto const normal{plane_.getNormal().getComponents()};
node["plane"]["normal"].push_back(normal.getX().magnitude());
node["plane"]["normal"].push_back(normal.getY().magnitude());
node["plane"]["normal"].push_back(normal.getZ().magnitude());
// the normal vector of the plane
auto const normal{plane_.getNormal().getComponents()};
node["plane"]["normal"].push_back(normal.getX().magnitude());
node["plane"]["normal"].push_back(normal.getY().magnitude());
node["plane"]["normal"].push_back(normal.getZ().magnitude());
// the x-axis vector
auto const xAxis_coords{xAxis_.getComponents(xAxis_.getCoordinateSystem())};
node["x-axis"].push_back(xAxis_coords.getX().magnitude());
node["x-axis"].push_back(xAxis_coords.getY().magnitude());
node["x-axis"].push_back(xAxis_coords.getZ().magnitude());
// the x-axis vector
auto const xAxis_coords{xAxis_.getComponents(xAxis_.getCoordinateSystem())};
node["x-axis"].push_back(xAxis_coords.getX().magnitude());
node["x-axis"].push_back(xAxis_coords.getY().magnitude());
node["x-axis"].push_back(xAxis_coords.getZ().magnitude());
// the y-axis vector
auto const yAxis_coords{yAxis_.getComponents(yAxis_.getCoordinateSystem())};
node["y-axis"].push_back(yAxis_coords.getX().magnitude());
node["y-axis"].push_back(yAxis_coords.getY().magnitude());
node["y-axis"].push_back(yAxis_coords.getZ().magnitude());
// the y-axis vector
auto const yAxis_coords{yAxis_.getComponents(yAxis_.getCoordinateSystem())};
node["y-axis"].push_back(yAxis_coords.getX().magnitude());
node["y-axis"].push_back(yAxis_coords.getY().magnitude());
node["y-axis"].push_back(yAxis_coords.getZ().magnitude());
node["delete_on_hit"] = deleteOnHit_;
node["delete_on_hit"] = deleteOnHit_;
return node;
}
return node;
}
template <typename TOutput>
void ObservationPlane<TOutput>::reset() {
......
......@@ -11,7 +11,7 @@
namespace corsika {
ObservationPlaneWriterParquet::ObservationPlaneWriterParquet()
: ParquetStreamer(){};
: ParquetStreamer(){}
void ObservationPlaneWriterParquet::startOfLibrary(
std::filesystem::path const& directory) {
......@@ -47,8 +47,8 @@ namespace corsika {
using namespace units::si;
// write the next row - we must write `shower_` first.
(*writer_) _ << shower_ << static_cast<int>(get_PDG(pid)) << energy / 1_eV << x / 1_m
<< y / 1_m << radius / 1_m << parquet::EndRow;
(*writer_) << shower_ << static_cast<int>(get_PDG(pid)) << energy / 1_eV << x / 1_m
<< y / 1_m << radius / 1_m << parquet::EndRow;
}
} // namespace corsika
......@@ -10,6 +10,6 @@
namespace corsika {
BaseOutput::BaseOutput()
: event_(0){}
: shower_(0){}
}
......@@ -11,8 +11,29 @@
#include <fstream>
#include <functional>
#include <iostream>
#include <iomanip>
#include <ctime>
#include <sstream>
namespace corsika {
std::string OutputManager::getCurrentTime() const {
// the format for our date string
auto fmt{"%d/%m/%Y %H:%M:%S %Z"};
// get the current time
auto t = std::time(nullptr);
auto current = *std::localtime(&t);
// create the string and push the time onto it
std::ostringstream oss;
oss << std::put_time(&current, fmt);
return oss.str();
}
void OutputManager::writeNode(YAML::Node const& node,
std::filesystem::path const& path) const {
......@@ -22,14 +43,11 @@ namespace corsika {
// and write the node to the output
out << node;
// open the output file - this is <output name>.json
// open the output file - this is <output name>.yaml
std::ofstream file(path.string());
// dump the JSON to the file
// dump the YAML to the file
file << out.c_str() << std::endl;
// and close the efile
file.close();
}
void OutputManager::writeTopLevelConfig() const {
......@@ -40,12 +58,25 @@ namespace corsika {
config["name"] = name_; // the simulation name
config["creator"] = "CORSIKA8"; // a tag to identify C8 libraries
config["version"] = "8.0.0-prealpha"; // the current version
// TODO: Add current datetime as a string to the config.
config["start time"] = getCurrentTime();
// write the node to a file
writeNode(config, root_ / ("config.yaml"));
}
// void OutputManager::writeTopLevelSummary() const {
// YAML::Node config;
// // some basic info
// config["start time"] = getCurrentTime(); // TODO:
// config["end time"] = getCurrentTime();
// // config["showers"] = 0; // TODO
// // write the node to a file
// writeNode(config, root_ / ("summary.yaml"));
// }
void OutputManager::initOutput(std::string const& name) const {
// construct the path to this directory
auto const path{root_ / name};
......@@ -77,7 +108,7 @@ namespace corsika {
throw std::runtime_error("Output directory already exists.");
}
// construct the directory for this run
// construct the directory for this library
std::filesystem::create_directory(root_);
// write the top level config file
......@@ -86,11 +117,19 @@ namespace corsika {
OutputManager::~OutputManager() {
// if we are being destructed but EndOfRun() has not been called,
// make sure that we gracefully close all the outputs
if (state_ == OutputState::EventInProgress || state_ == OutputState::RunInitialized) {
endOfRun();
if (state_ == OutputState::ShowerInProgress) {
// if this the destructor is called before the shower has been explicitly
// ended, print a warning and end the shower before continuing.
logger->warn(
"OutputManager was destroyed before endOfShower() called."
" The last shower in this libray may be incomplete.");
endOfShower();
}
// if we are being destructed but EndOfLibrary() has not been called,
// make sure that we gracefully close all the outputs. This is a supported
// method of operation so we don't issue a warning here
if (state_ == OutputState::LibraryReady) { endOfLibrary(); }
}
// void OutputManager::add(std::string const& name, BaseOutput& output) {
......@@ -112,15 +151,14 @@ namespace corsika {
initOutput(name);
}
void OutputManager::startOfRun() {
void OutputManager::startOfLibrary() {
// this is only valid when we haven't started a run
// or have already finished a run
if
(!(state_ == OutputState::RunNoInit || state_ == OutputState::RunFinished)) {
// this is only valid when we haven't started a library
// or have already finished a library
if (!(state_ == OutputState::NoInit || state_ == OutputState::LibraryFinished)) {
throw std::runtime_error("startOfRun() called in invalid state.");
}
throw std::runtime_error("startOfLibrary() called in invalid state.");
}
// we now forward this signal to all of our outputs
for (auto& [name, output] : outputs_) {
......@@ -128,49 +166,47 @@ namespace corsika {
// construct the path to this output subdirectory
auto const path{root_ / name};
// and start the run
output.get().startOfRun(path);
// and start the library
output.get().startOfLibrary(path);
}
// we have now started running
state_ = OutputState::RunInitialized;
state_ = OutputState::LibraryReady;
}
void OutputManager::startOfEvent() {
void OutputManager::startOfShower() {
// if this is called but we are still in the initialized state,
// make sure that we transition to EventInProgress
// if (state_ == OutputState::RunNoInit) { startOfRun(); }
// if this is called and we still in the "no init" state, then
// this is the first shower in the library so make sure we start it
if (state_ == OutputState::NoInit) { startOfLibrary(); }
// now start the event for all the outputs
for (auto& [name, output] : outputs_) { output.get().startOfEvent(); }
for (auto& [name, output] : outputs_) { output.get().startOfShower(); }
// and transition to the in progress state
state_ = OutputState::EventInProgress;
state_ = OutputState::ShowerInProgress;
}
void OutputManager::endOfEvent() {
void OutputManager::endOfShower() {
for (auto& [name, output] : outputs_) { output.get().endOfEvent(); }
for (auto& [name, output] : outputs_) { output.get().endOfShower(); }
// switch back to the initialized state
state_ = OutputState::RunInitialized;
state_ = OutputState::LibraryReady;
}
void OutputManager::endOfRun() {
void OutputManager::endOfLibrary() {
// we can only call endOfRun when we have already started
if (state_ == OutputState::RunNoInit) {
throw std::runtime_error("endOfRun() called in invalid state.");
// we can only call endOfLibrary when we have already started
if (state_ == OutputState::NoInit) {
throw std::runtime_error("endOfLibrary() called in invalid state.");
}
for (auto& [name, output] : outputs_) { output.get().endOfRun(); }
// and the run has finished
state_ = OutputState::RunFinished;
// forward the endOfLibrary() call to all the registered outputs
for (auto& [name, output] : outputs_) { output.get().endOfLibrary(); }
// write any final state information into the config files
// for (auto& [name, output] : outputs_) { output.get().endOfRun(); }
// and the library has finished
state_ = OutputState::LibraryFinished;
}
} // namespace corsika
......@@ -21,9 +21,7 @@ namespace corsika {
builder_.created_by("CORSIKA8");
// add run and event tags to the file
addField("run", parquet::Repetition::REQUIRED, parquet::Type::INT32,
parquet::ConvertedType::INT_32);
addField("event", parquet::Repetition::REQUIRED, parquet::Type::INT32,
addField("shower", parquet::Repetition::REQUIRED, parquet::Type::INT32,
parquet::ConvertedType::INT_32);
}
......
......@@ -55,4 +55,4 @@ namespace corsika {
} // namespace corsika
#include <corsika/details/modules/writers/ObservationPlaneWriterParquet.inl>
#include <corsika/detail/modules/writers/ObservationPlaneWriterParquet.inl>
......@@ -18,12 +18,10 @@ namespace corsika {
* This is the base class for all outputs so that they
* can be stored in homogeneous containers.
*/
// class BaseOutput : public std::enable_shared_from_this<BaseOutput> {
class BaseOutput {
protected:
int event_{0}; ///< The current event number.
int run_{0}; ///< The current run number.
int shower_{0}; ///< The current event number.
BaseOutput();
......@@ -31,22 +29,22 @@ namespace corsika {
/**
* Called at the start of each run.
*/
virtual void startOfRun(std::filesystem::path const& directory) = 0;
virtual void startOfLibrary(std::filesystem::path const& directory) = 0;
/**
* Called at the start of each event/shower.
*/
virtual void startOfEvent() {}
virtual void startOfShower() {}
/**
* Called at the end of each event/shower.
*/
virtual void endOfEvent() = 0;
virtual void endOfShower() = 0;
/**
* Called at the end of each run.
*/
virtual void endOfRun() = 0;
virtual void endOfLibrary() = 0;
/**
* Get the configuration of this output.
......
......@@ -23,21 +23,20 @@ namespace corsika {
* Indicates the current state of this manager.
*/
enum class OutputState {
RunNoInit,
RunInitialized,
EventInProgress,
RunFinished,
NoInit,
LibraryReady,
ShowerInProgress,
LibraryFinished,
};
OutputState state_{OutputState::RunNoInit}; ///< The current state of this manager.
std::string const name_; ///< The name of this simulation file.
std::filesystem::path const root_; ///< The top-level directory for the output.
OutputState state_{OutputState::NoInit}; ///< The current state of this manager.
std::string const name_; ///< The name of this simulation file.
std::filesystem::path const root_; ///< The top-level directory for the output.
inline static auto logger{get_logger("output")}; ///< A custom logger.
/**
* The outputs that have been registered with us.
*/
// std::map<std::string, std::shared_ptr<BaseOutput>> outputs_;
std::map<std::string, std::reference_wrapper<BaseOutput>> outputs_;
/**
......@@ -55,6 +54,11 @@ namespace corsika {
*/
void initOutput(std::string const& name) const;
/**
* Get the current local time as a string.
*/
std::string getCurrentTime() const;
public:
/**
* Construct an OutputManager instance with a name in a given directory.
......@@ -78,32 +82,31 @@ namespace corsika {
*/
template <typename TOutput>
void add(std::string const& name, TOutput& output);
// void add(std::string const& name, BaseOutput& output);
/**
* Called at the start of each run.
* Called at the start of each library.
*
* This iteratively calls startOfRun on each registered output.
* This iteratively calls startOfLibrary on each registered output.
*/
void startOfRun();
void startOfLibrary();
/**
* Called at the start of each event/shower.
* This iteratively calls startOfEvent on each registered output.
*/
void startOfEvent();
void startOfShower();
/**
* Called at the end of each event/shower.
* This iteratively calls endOfEvent on each registered output.
*/
void endOfEvent();
void endOfShower();
/**
* Called at the end of each run.
* This iteratively calls endOfRun on each registered output.
* Called at the end of each library.
* This iteratively calls endOfLibrary on each registered output.
*/
void endOfRun();
void endOfLibrary();
}; // class OutputManager
......
......@@ -190,9 +190,6 @@ int main(int argc, char** argv) {
ShowerAxis const showerAxis{injectionPos, (showerCore - injectionPos) * 1.5, env};
// create the output manager that we then register outputs with
OutputManager output("vertical_EAS_outputs");
// setup processes, decays and interactions
corsika::sibyll::Interaction sibyll;
......@@ -239,6 +236,9 @@ int main(int argc, char** argv) {
Plane const obsPlane(showerCore, DirectionVector(rootCS, {0., 0., 1.}));
ObservationPlane observationLevel(obsPlane, DirectionVector(rootCS, {1., 0., 0.}));
// create the output manager that we then register outputs with
OutputManager output("vertical_EAS_outputs");
// register the observation plane with the output
output.add("obsplane", observationLevel);
......@@ -272,9 +272,7 @@ int main(int argc, char** argv) {
// to fix the point of first interaction, uncomment the following two lines:
// EAS.forceInteraction();
output.startOfRun();
EAS.run();
output.endOfRun();
cut.showResults();
em_continuous.showResults();
......@@ -294,4 +292,6 @@ int main(int argc, char** argv) {
save_hist(hists.labHist(), "inthist_lab_verticalEAS.npz", true);
save_hist(hists.CMSHist(), "inthist_cms_verticalEAS.npz", true);
longprof.save("longprof_verticalEAS.txt");
output.endOfLibrary();
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment