diff --git a/Documentation/Examples/vertical_EAS.cc b/Documentation/Examples/vertical_EAS.cc index 388255bfc6d5dc286173fa3fb31ffc1a2115963e..0f4d57c8787dfcdf444409075732c28457fe3140 100644 --- a/Documentation/Examples/vertical_EAS.cc +++ b/Documentation/Examples/vertical_EAS.cc @@ -162,8 +162,20 @@ int main(int argc, char** argv) { stack.AddParticle(std::make_tuple(beamCode, E0, plab, injectionPos, 0_ns, A, Z)); } else { - stack.AddParticle( - std::make_tuple(particles::Code::Proton, E0, plab, injectionPos, 0_ns)); + if (Z == 1) { + stack.AddParticle(std::tuple<particles::Code, units::si::HEPEnergyType, + corsika::stack::MomentumVector, geometry::Point, + units::si::TimeType>{particles::Code::Proton, E0, plab, + injectionPos, 0_ns}); + } else if (Z == 0) { + stack.AddParticle(std::tuple<particles::Code, units::si::HEPEnergyType, + corsika::stack::MomentumVector, geometry::Point, + units::si::TimeType>{particles::Code::Neutron, E0, + plab, injectionPos, 0_ns}); + } else { + std::cerr << "illegal parameters" << std::endl; + return EXIT_FAILURE; + } } // we make the axis much longer than the inj-core distance since the @@ -217,8 +229,8 @@ int main(int argc, char** argv) { process::longitudinal_profile::LongitudinalProfile longprof{showerAxis}; Plane const obsPlane(showerCore, Vector<dimensionless_d>(rootCS, {0., 0., 1.})); - process::observation_plane::ObservationPlane observationLevel(obsPlane, - "particles.dat"); + process::observation_plane::ObservationPlane observationLevel( + obsPlane, Vector<dimensionless_d>(rootCS, {1., 0., 0.}), "particles.dat"); process::UrQMD::UrQMD urqmd; process::interaction_counter::InteractionCounter urqmdCounted{urqmd}; diff --git a/Processes/ObservationPlane/ObservationPlane.cc b/Processes/ObservationPlane/ObservationPlane.cc index ab3a8023215528084dd242b656597f74e9989d23..7abff35139b375ccdf7e82fb50287730695351d9 100644 --- a/Processes/ObservationPlane/ObservationPlane.cc +++ b/Processes/ObservationPlane/ObservationPlane.cc @@ -14,14 +14,18 @@ using namespace corsika::process::observation_plane; using namespace corsika::units::si; -ObservationPlane::ObservationPlane(geometry::Plane const& obsPlane, - std::string const& filename, bool deleteOnHit) +ObservationPlane::ObservationPlane( + geometry::Plane const& obsPlane, + geometry::Vector<units::si::dimensionless_d> const& x_axis, + std::string const& filename, bool deleteOnHit) : plane_(obsPlane) , outputStream_(filename) , deleteOnHit_(deleteOnHit) , energy_ground_(0_GeV) , count_ground_(0) { - outputStream_ << "#PDG code, energy / eV, distance to center / m" << std::endl; + , xAxis_(x_axis.normalized()) + , yAxis_(obsPlane.GetNormal().cross(xAxis_)) { + outputStream_ << "#PDG code, energy / eV, x distance / m, y distance / m" << std::endl; } corsika::process::EProcessReturn ObservationPlane::DoContinuous( @@ -37,8 +41,11 @@ corsika::process::EProcessReturn ObservationPlane::DoContinuous( } const auto energy = particle.GetEnergy(); + auto const displacement = trajectory.GetPosition(1) - plane_.GetCenter(); + outputStream_ << static_cast<int>(particles::GetPDG(particle.GetPID())) << ' ' << energy / 1_eV << ' ' + << displacement.dot(xAxis_) / 1_m << ' ' << displacement.dot(yAxis_) / 1_m << (trajectory.GetPosition(1) - plane_.GetCenter()).norm() / 1_m << std::endl; diff --git a/Processes/ObservationPlane/ObservationPlane.h b/Processes/ObservationPlane/ObservationPlane.h index 86c5ba0f2100043f5a7c4b9c2ba250b1e3a57a20..3b530f6213a6d88160329c3165e004a3ec6b66b7 100644 --- a/Processes/ObservationPlane/ObservationPlane.h +++ b/Processes/ObservationPlane/ObservationPlane.h @@ -26,7 +26,9 @@ namespace corsika::process::observation_plane { class ObservationPlane : public corsika::process::ContinuousProcess<ObservationPlane> { public: - ObservationPlane(geometry::Plane const&, std::string const&, bool = true); + ObservationPlane(geometry::Plane const&, + geometry::Vector<units::si::dimensionless_d> const&, + std::string const&, bool = true); corsika::process::EProcessReturn DoContinuous( corsika::setup::Stack::ParticleType& vParticle, @@ -47,5 +49,6 @@ namespace corsika::process::observation_plane { units::si::HEPEnergyType energy_ground_ = 0 * units::si::electronvolt; unsigned int count_ground_ = 0; + geometry::Vector<units::si::dimensionless_d> const xAxis_, yAxis_; }; } // namespace corsika::process::observation_plane