diff --git a/corsika/detail/framework/geometry/Point.inl b/corsika/detail/framework/geometry/Point.inl
index 156723dd149102ecec8ee88ab773dafa5cfec84d..97f932e2ddb942af8e9c88a67bd7bd86cbe35d90 100644
--- a/corsika/detail/framework/geometry/Point.inl
+++ b/corsika/detail/framework/geometry/Point.inl
@@ -35,6 +35,10 @@ namespace corsika {
return getCoordinates(pCS).getZ();
}
+ inline LengthType Point::distance_to(Point const& point) const {
+ return (*this - point).getNorm();
+ }
+
/// this always returns a QuantityVector as triple
inline QuantityVector<length_d> Point::getCoordinates(
CoordinateSystemPtr const& pCS) const {
diff --git a/corsika/modules/radio/CoREAS.hpp b/corsika/modules/radio/CoREAS.hpp
new file mode 100755
index 0000000000000000000000000000000000000000..ca1ea7e6266e3acd03dbd34a55d628f0a9f32445
--- /dev/null
+++ b/corsika/modules/radio/CoREAS.hpp
@@ -0,0 +1,67 @@
+/*
+ * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+#pragma once
+
+#include <corsika/modules/radio/RadioProcess.hpp>
+#include <corsika/modules/radio/propagators/StraightPropagator.hpp>
+#include <corsika/framework/geometry/QuantityVector.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+
+namespace corsika {
+
+ /**
+ *
+ */
+ template <typename TPropagator = StraightPropagator>
+ class CoREAS final : public RadioProcess<CoREAS> {
+
+
+ public:
+ /**
+ * Simulate the radio emission from a particle across a track.
+ *
+ * This must be provided by the TRadioImpl.
+ *
+ * @param particle The current particle.
+ * @param track The current track.
+ *
+ */
+ template <typename Particle, typename Track>
+ EProcessReturn simulate(Particle&, Track const&) const {
+
+ // loop over every antenna
+ for (auto& antenna : getAntennas()) {
+
+ // use TPropagator to calculate Path from track to antenna
+ auto paths = ;/* a collection of paths */
+
+ // do endpoint + ZHS formalism
+
+ // give the electric field and the direction to the antenna
+ antenna.receive(...); // << something
+ }
+
+ }
+
+ /**
+ * Return the maximum step length for this particle and track.
+ *
+ * This must be provided by the TRadioImpl.
+ *
+ * @param particle The current particle.
+ * @param track The current track.
+ *
+ * @returns The maximum length of this track.
+ */
+ template <typename Particle, typename Track>
+ LengthType MaxStepLength(Particle const& particle,
+ Track const& track) const;
+
+ }; // END: class RadioProcess
+
+} // namespace corsika
diff --git a/corsika/modules/radio/RadioProcess.hpp b/corsika/modules/radio/RadioProcess.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..03f93f9f50e2132589b9aa1f28a01a9de7ea9d83
--- /dev/null
+++ b/corsika/modules/radio/RadioProcess.hpp
@@ -0,0 +1,84 @@
+/*
+ * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+#pragma once
+
+#include <corsika/framework/process/ContinuousProcess.hpp>
+
+namespace corsika {
+
+ /**
+ * The base interface for radio emission processes.
+ *
+ * TRadioImpl is the concrete implementation of the radio algorithm.
+ * TRadioDetector is the detector instance that stores antennas
+ * and is responsible for managing the output writing.
+ */
+ template <typename TRadioDetector, typename TRadioImpl, typename TPropagator>
+ class RadioProcess : public ContinuousProcess<
+ RadioProcess<TRadioDetector, TRadioImpl, TPropagator>> {
+
+ /**
+ * Get a reference to the underlying radio implementation.
+ */
+ TRadioImpl& implementation() { return static_cast<TRadioImpl&>(*this); }
+
+ /*
+ * Get a const reference to the underlying implementation.
+ */
+ TRadioImpl const& implementation() const {
+ return static_cast<TRadioImpl const&>(*this);
+ }
+
+ protected:
+ TRadioDetector& detector_; ///< The radio detector we store into.
+ TPropagator propagator_; ///< The propagator implementation.
+
+ public:
+ /**
+ * Construct a new RadioProcess.
+ */
+ template <typename... TArgs>
+ RadioProcess(TRadioDetector& detector, TArgs&&... args)
+ : detector_(detector)
+ , propagator_(args...) {}
+
+ /**
+ * Perform the continuous process (radio emission).
+ *
+ * This handles filtering individual particle tracks
+ * before passing them to `Simulate`.`
+ *
+ * @param particle The current particle.
+ * @param track The current track.
+ */
+ template <typename Particle, typename Track>
+ ProcessReturn DoContinuous(Particle& particle, Track const& track) const {
+ // we wrap Simulate() in DoContinuous as the plan is to add particle level
+ // filtering or thinning for calculation of the radio emission. This is
+ // important for controlling the runtime of radio (by ignoring particles
+ // that aren't going to contribute i.e. heavy hadrons)
+ return this->implementation().simulate(particle, track);
+ }
+
+ /**
+ * TODO: This is placeholder so we can use text output while
+ * we wait for the true output formatting to be ready.
+ **/
+ bool writeOutput() const {
+
+ // how this method should work:
+ // 1. Loop over the antennas in the collection
+ // 2. Get their waveforms
+ // 3. Create a textfile for each antenna
+ // 4. and write out two columns, time and field.
+
+ }
+
+ }; // END: class RadioProcess
+
+} // namespace corsika::process::radio
diff --git a/corsika/modules/radio/antennas/Antenna.hpp b/corsika/modules/radio/antennas/Antenna.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..8596d576334b6a8f66bddd6a0c5ad647b5d570f3
--- /dev/null
+++ b/corsika/modules/radio/antennas/Antenna.hpp
@@ -0,0 +1,79 @@
+/*
+ * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+#pragma once
+
+#include <corsika/framework/geometry/Point.hpp>
+
+namespace corsika {
+
+ /**
+ * A common abstract interface for radio antennas.
+ *
+ * All concrete antenna implementations should be of
+ * type Antenna<T> where T is a concrete antenna implementation.
+ *
+ */
+ template <typename AntennaImpl>
+ class Antenna {
+
+ std::string const name_; ///< The name/identifier of this antenna.
+ Point const location_; ///< The location of this antenna.
+
+ public:
+ // this stores the polarization vector of an electric field
+ using ElectricFieldVector =
+ QuantityVector<ElectricFieldType>;
+ using MagneticFieldVector =
+ QuantityVector<MagneticFieldType>;
+
+ // a dimensionless vector used for the incident direction
+ using Vector = QuantityVector<dimensionless_d>;
+
+ /**
+ * \brief Construct a base antenna instance.
+ *
+ * @param name A name for this antenna.
+ * @param location The location of this antenna.
+ *
+ */
+ Antenna(std::string const& name, Point const& location)
+ : name_(name)
+ , location_(location){};
+
+ /**
+ * Receive a signal at this antenna.
+ *
+ * This is a general implementation call that must be specialized
+ * for the particular antenna implementation and usage.
+ *
+ */
+ // template <typename... TVArgs>
+ // void Receive(TVArgs&& args...);
+
+ /**
+ * Get the location of this antenna.
+ */
+ Point const& getLocation() const { return location_; };
+
+ /**
+ * Get the name of this name antenna.
+ *
+ * This is used in producing the output data file.
+ */
+ std::string const& getName() const { return name_; };
+
+ /**
+ * Reset the antenna before starting a new simulation.
+ */
+ void reset();
+
+ }; // END: class Antenna final
+
+} // namespace corsika
diff --git a/corsika/modules/radio/antennas/TimeDomainAntenna.hpp b/corsika/modules/radio/antennas/TimeDomainAntenna.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..4ce0c681e6acb8934aae5d55038c08b104e4eff9
--- /dev/null
+++ b/corsika/modules/radio/antennas/TimeDomainAntenna.hpp
@@ -0,0 +1,86 @@
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+#pragma once
+
+#include <corsika/modules/radio/antennas/Antenna.hpp>
+
+namespace corsika {
+
+ /**
+ * An implementation of a time-domain antenna that has a customized
+ * start time, sampling period, and waveform duration.
+ *
+ */
+ class TimeDomainAntenna final : public Antenna<TimeDomainAntenna> {
+
+ using Array = std::vector<double>; ///< The type used to store arrays.
+
+ TimeType const start_time_; ///< The start time of this waveform.
+ TimeType const duration_; ///< The duration of this waveform.
+ TimeType const sampling_period_; ///< The sampling period of this antenna.
+
+ Array waveform_; ///< The waveform stored by this antenna.
+
+ protected:
+ // expose the CRTP interfaces constructor
+
+ public:
+ // import the methods from the antenna
+ using Antenna<TimeDomainAntenna>::Antenna;
+ using Antenna<TimeDomainAntenna>::getName;
+
+ /**
+ * Construct a new TimeDomainAntenna.
+ *
+ * @param name The name of this antenna.
+ * @param location The location of this antenna.
+ * @param start_time The starting time of this waveform.
+ * @param duration The duration of this waveform.
+ * @param sampling_period The sampling period of this waveform.
+ *
+ */
+ TimeDomainAntenna(std::string const& name, Point const& location,
+ TimeType const& start_time,
+ TimeType const& duration,
+ TimeType const& sampling_period)
+ : Antenna(name, location)
+ , start_time_(start_time)
+ , duration_(duration)
+ , sampling_period_(sampling_period){};
+
+ /**
+ * Receive an electric field at this antenna.
+ *
+ * This assumes that the antenna will receive
+ * an *instantaneous* electric field modeled as a delta function.
+ *
+ * @param time The (global) time at which this
+ * @param field The incident electric field vector.
+ *
+ */
+ void receive(TimeType const time, ElectricFieldVector const& efield) const;
+
+ /**
+ * Get the current waveform for this antenna.
+ *
+ * NOTE: Currently returns ns and V/m but this should be UNITful
+ *
+ * @returns A pair of the sample times, and the field
+ */
+ std::pair<Array, Array> getWaveform() const;
+
+ /**
+ * Reset the antenna before starting a new simulation.
+ */
+ void reset() { waveform_.clear(); };
+
+ }; // END: class Antenna final
+
+} // namespace corsika
diff --git a/corsika/modules/radio/detectors/RadioDetector.hpp b/corsika/modules/radio/detectors/RadioDetector.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..7c21d7e48631a1f63b733177988b136d47db9ffa
--- /dev/null
+++ b/corsika/modules/radio/detectors/RadioDetector.hpp
@@ -0,0 +1,50 @@
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+#pragma once
+
+namespace corsika {
+
+ /**
+ * The base interface for radio detectors.
+ */
+
+ template <typename TAntennaImpl, typename TDetectorImpl>
+ class RadioDetector {
+
+ /**
+ * The collection of antennas used in this simulation.
+ */
+ std::vector<TAntennaImpl> antennas_;
+
+ public:
+ /**
+ * Add an antenna to this radio process.
+ *
+ * @param antenna The antenna to add
+ */
+ auto addAntenna(TAntennaImpl const& antenna) -> void { antennas_.push_back(antenna); }
+
+ /**
+ * Get a *non*-const reference to the collection of antennas.
+ *
+ * @returns An iterable mutable reference to the antennas.
+ */
+ std::vector<TAntennaImpl> const& getAntennas() { return antennas_; }
+
+ /**
+ * Reset all the antenna waveforms.
+ */
+ auto reset() -> void {
+ std::for_each(antennas_.begin(), antennas_.end(), std::mem_fn(&TAntennaImpl::Reset));
+ };
+
+ }; // END: class RadioDetector
+
+} // namespace corsika
diff --git a/corsika/modules/radio/detectors/TimeDomainDetector.hpp b/corsika/modules/radio/detectors/TimeDomainDetector.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..23974262f92b63d5eb1417110aa18712ae295dc2
--- /dev/null
+++ b/corsika/modules/radio/detectors/TimeDomainDetector.hpp
@@ -0,0 +1,27 @@
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+#pragma once
+
+#include <vector>
+
+#include "corsika/modules/radio/detectors/RadioDetector.hpp"
+
+namespace corsika {
+
+ /**
+ * A common interface for radio detectors.
+ */
+ template <typename TAntennaImpl>
+ class TimeDomainDetector final : public RadioDetector<TAntennaImpl,TimeDomainDetector<TAntennaImpl>> {
+
+ public:
+ }; // END: class RadioDetector
+
+} // namespace corsika
diff --git a/corsika/modules/radio/propagators/RadioPropagator.hpp b/corsika/modules/radio/propagators/RadioPropagator.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..394e2570c30217c2887e5b2b8d44a52b7fc9e492
--- /dev/null
+++ b/corsika/modules/radio/propagators/RadioPropagator.hpp
@@ -0,0 +1,46 @@
+/*
+ * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+#pragma once
+
+#include <vector>
+
+#include <corsika/framework/geometry/Path.hpp>
+#include <corsika/modules/radio/propagators/SignalPath.hpp>
+
+namespace corsika {
+
+ /**
+ * Radio propagators are used to calculate the propagation
+ * paths from particles to antennas. Any class that wants
+ * to be used as a RadioPropagator must implement the
+ * following methods:
+ *
+ * SignalPathCollection Propagate(Point const& start,
+ * Point const& end,
+ * LengthType const stepsize);
+ */
+ template <typename TImpl, typename TEnvironment>
+ class RadioPropagator {
+
+ protected:
+ // Since we typically know roughly how many paths will
+ // be computed for a given propagator, we use a std::vector here.
+ using SignalPathCollection = std::vector<SignalPath> const;
+
+ TEnvironment const& env_; ///< The environment.
+
+ public:
+ /**
+ * Construct a new RadioPropagator instance.
+ */
+ RadioPropagator(TEnvironment const& env)
+ : env_(env) {}
+
+ }; // class RadioPropagator
+
+} // namespace corsika
diff --git a/corsika/modules/radio/propagators/SignalPath.hpp b/corsika/modules/radio/propagators/SignalPath.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..76157d899473b5fd34d53bc9c8e47a1299fe8e6b
--- /dev/null
+++ b/corsika/modules/radio/propagators/SignalPath.hpp
@@ -0,0 +1,45 @@
+/*
+ * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
+#pragma once
+
+#include <corsika/framework/geometry/Path.hpp>
+#include <corsika/framework/geometry/Vector.hpp>
+
+namespace corsika {
+
+ /**
+ * Store the photon signal path between two points.
+ *
+ * This is basically a container class
+ */
+ struct SignalPath final : private Path {
+
+ using path = std::deque<Point>;
+
+ TimeType const time_; ///< The total propagation time.
+ double const average_refractivity_; ///< The average refractivity.
+ Vector<dimensionless_d> const emit_; ///< The (unit-length) emission vector.
+ Vector<dimensionless_d> const receive_; ///< The (unit-length) receive vector.
+ path const points_; ///< A collection of points that make up the geometrical path.
+
+ /**
+ * Create a new SignalPath instance.
+ */
+ SignalPath(TimeType const time, double const average_refractivity,
+ Vector<dimensionless_d> const emit, Vector<dimensionless_d> const receive,
+ path const& points)
+ : Path(points)
+ , time_(time)
+ , average_refractivity_(average_refractivity)
+ , emit_(emit)
+ , receive_(receive) {}
+
+ }; // class SignalPath
+
+} // namespace corsika
diff --git a/corsika/modules/radio/propagators/StraightPropagator.hpp b/corsika/modules/radio/propagators/StraightPropagator.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..3ee2f3bf075c795766b30afacc82908de0293860
--- /dev/null
+++ b/corsika/modules/radio/propagators/StraightPropagator.hpp
@@ -0,0 +1,129 @@
+/*
+ * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+#pragma once
+
+#include <corsika/media/Environment.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/geometry/Vector.hpp>
+#include <corsika/framework/core/PhysicalConstants.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/modules/radio/propagators/RadioPropagator.hpp>
+
+namespace corsika {
+
+ /**
+ * This class implements a basic propagator that uses
+ * the straight-line (vector) between the particle
+ * location and the antenna as the trajectory.
+ *
+ * This is what is used in ZHAireS and CoREAS in C7.
+ */
+ template <typename TEnvironment>
+ class StraightPropagator final
+ : public RadioPropagator<StraightPropagator<TEnvironment>, TEnvironment> {
+
+ using Base = RadioPropagator<StraightPropagator<TEnvironment>, TEnvironment>;
+ using SignalPathCollection = typename Base::SignalPathCollection;
+
+ public:
+
+ /**
+ * Construct a new StraightPropagator with a given environment.
+ *
+ */
+ StraightPropagator(TEnvironment const& env)
+ : RadioPropagator<StraightPropagator, TEnvironment>(env){};
+ // TODO: maybe the constructor doesn't take any arguments for the environment (?)
+
+ /**
+ * Return the collection of paths from `start` to `end`.
+ * or from 'source' which is the emission point to 'destination'
+ * which is the location of the antenna
+ */
+ SignalPathCollection propagate(Point const& source,
+ Point const& destination,
+ LengthType const stepsize) const {
+
+ /**
+ * get the normalized (unit) vector from `source` to `destination'.
+ * this is also the `emit` and `receive` vectors in the SignalPath class.
+ * in this case emit and receive unit vectors should be the same
+ * so they are both called direction
+ */
+ auto direction{(destination - source).normalized()};
+
+ // the step is the direction vector with length `stepsize`
+ auto step{direction * stepsize};
+
+ //calculate the number of points (roughly) for the numerical integration.
+ auto n_points {(destination - source).getNorm() / stepsize};
+
+ // get the universe for this environment
+ auto const* const universe{Base::env_.getUniverse().get()};
+
+ // the points that we build along the way for the numerical integration
+ std::deque<Point> points;
+
+ // store value of the refractive index at points
+ std::vector<double> rindex;
+ rindex.reserve(n_points);
+
+ // loop from `source` to `destination` to store values before Simpson's rule.
+ // this loop skips the last point 'destination'
+ for (auto point = source; (point - destination).getNorm() > 0.6 * stepsize;
+ point = point + step) {
+
+ // get the environment node at this specific 'point'
+ auto const* node{universe->getContainingNode(point)};
+
+ // get the associated refractivity at 'point'
+ auto const refractive_index{node->getModelProperties().getRefractiveIndex(point)};
+ rindex.push_back(refractive_index);
+
+ // add this 'point' to our deque collection
+ points.push_back(point);
+ }
+
+ //add the refractive index of last point 'destination' and store it
+ auto const* node{universe->getContainingNode(destination)};
+ auto const refractive_index{node->getModelProperties().getRefractiveIndex(destination)};
+ rindex.push_back(refractive_index);
+ points.push_back(destination);
+
+ // Apply Simpson's rule
+ auto N = rindex.size();
+ std::size_t index = 0;
+ double sum = rindex.at(index);
+ auto refra_ = rindex.at(index);
+ auto h = ((destination - source).getNorm()) / (N - 1);
+ for (std::size_t index = 1; index < (N - 1); index += 2) {
+ sum += 4 * rindex.at(index);
+ refra_ += rindex.at(index);
+ }
+ for (std::size_t index = 2; index < (N - 1); index += 2) {
+ sum += 2 * rindex.at(index);
+ refra_ += rindex.at(index);
+ }
+ index = N - 1;
+ sum = sum + rindex.at(index);
+ refra_ += rindex.at(index);
+
+ // compute the total time delay.
+ TimeType time = sum * (h / (3 * constants::c));
+
+ // compute the average refractivity.
+ auto average_refractivity = refra_ / N;
+
+ // realize that emission and receive vector are 'direction' in this case.
+ return { SignalPath(time, average_refractivity, direction , direction, points) };
+
+ } // END: propagate()
+
+ }; // End: StraightPropagator
+
+} // namespace corsika
\ No newline at end of file
diff --git a/tests/modules/testRadio.cpp b/tests/modules/testRadio.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..7b7321e86ad986453885afa2bf403baeb4c3a9dd
--- /dev/null
+++ b/tests/modules/testRadio.cpp
@@ -0,0 +1,280 @@
+/*
+ * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+#include <catch2/catch.hpp>
+
+#include <corsika/modules/radio/ZHS.hpp>
+#include <corsika/modules/radio/antennas/TimeDomainAntenna.hpp>
+#include <corsika/modules/radio/detectors/TimeDomainDetector.hpp>
+#include <corsika/modules/radio/propagators/StraightPropagator.hpp>
+#include <corsika/modules/radio/propagators/SignalPath.hpp>
+#include <corsika/modules/radio/propagators/RadioPropagator.hpp>
+
+
+#include <corsika/media/Environment.hpp>
+#include <corsika/media/HomogeneousMedium.hpp>
+#include <corsika/media/IMediumModel.hpp>
+#include <corsika/media/IRefractiveIndexModel.hpp>
+#include <corsika/media/LayeredSphericalAtmosphereBuilder.hpp>
+#include <corsika/media/UniformRefractiveIndex.hpp>
+#include <corsika/media/ExponentialRefractiveIndex.hpp>
+#include <corsika/media/VolumeTreeNode.hpp>
+#include <corsika/framework/geometry/CoordinateSystem.hpp>
+#include <corsika/framework/geometry/Line.hpp>
+#include <corsika/framework/geometry/Point.hpp>
+#include <corsika/framework/geometry/RootCoordinateSystem.hpp>
+#include <corsika/framework/geometry/Vector.hpp>
+#include <corsika/setup/SetupStack.hpp>
+#include <corsika/setup/SetupTrajectory.hpp>
+#include <corsika/framework/core/PhysicalUnits.hpp>
+#include <corsika/framework/core/PhysicalConstants.hpp>
+
+using namespace corsika;
+
+double constexpr absMargin = 1.0e-7;
+
+TEST_CASE("Radio", "[processes]") {
+
+// logging::set_level(logging::level::info);
+// corsika_logger->set_pattern("[%n:%^%-8l%$] custom pattern: %v");
+
+
+ // check that I can create and reset a TimeDomain process
+ SECTION("TimeDomainDetector") {
+
+ // construct a time domain detector
+ TimeDomainDetector<TimeDomainAntenna> detector;
+
+ // // and construct some time domain radio process
+ // TimeDomain<ZHS<>, TimeDomainDetector<TimeDomainAntenna>> process(detector);
+ // TimeDomain<ZHS<CPU>, TimeDomainDetector<TimeDomainAntenna>> process2(detector);
+ // TimeDomain<ZHS<CPU>, decltype(detector)> process3(detector);
+ }
+
+ // check that I can create time domain antennas
+ SECTION("TimeDomainDetector") {
+
+ // create an environment so we can get a coordinate system
+ Environment<IMediumModel> env;
+
+ // the antenna location
+ const auto point{Point(env.getCoordinateSystem(), 1_m, 2_m, 3_m)};
+
+ // check that I can create an antenna at (1, 2, 3)
+ const auto ant{TimeDomainAntenna("antenna_name", point)};
+
+ // assert that the antenna name is correct
+ REQUIRE(ant.getName() == "antenna_name");
+
+ // and check that the antenna is at the right location
+ REQUIRE((ant.getLocation() - point).getNorm() < 1e-12 * 1_m);
+
+ // construct a radio detector instance to store our antennas
+ TimeDomainDetector<TimeDomainAntenna> detector;
+
+ // add this antenna to the process
+ detector.addAntenna(ant);
+ }
+
+ // check that I can create working Straight Propagators in different environments
+ SECTION("Straight Propagator w/ Uniform Refractive Index") {
+
+ // create an environment with uniform refractive index of 1
+ using UniRIndex =
+ UniformRefractiveIndex<HomogeneousMedium<IRefractiveIndexModel<IMediumModel>>>;
+
+ using EnvType = Environment<IRefractiveIndexModel<IMediumModel>>;
+ EnvType env;
+
+ // get a coordinate system
+ const CoordinateSystemPtr rootCS = env.getCoordinateSystem();
+
+ auto Medium = EnvType::createNode<Sphere>(
+ Point{rootCS, 0_m, 0_m, 0_m}, 1_km * std::numeric_limits<double>::infinity());
+
+ auto const props = Medium->setModelProperties<UniRIndex>(
+ 1, 1_kg / (1_m * 1_m * 1_m),
+ NuclearComposition(
+ std::vector<Code>{Code::Nitrogen},
+ std::vector<float>{1.f}));
+
+ env.getUniverse()->addChild(std::move(Medium));
+
+ // get some points
+ Point p0(rootCS, {0_m, 0_m, 0_m});
+ // Point p1(rootCS, {0_m, 0_m, 1_m});
+ // Point p2(rootCS, {0_m, 0_m, 2_m});
+ // Point p3(rootCS, {0_m, 0_m, 3_m});
+ // Point p4(rootCS, {0_m, 0_m, 4_m});
+ // Point p5(rootCS, {0_m, 0_m, 5_m});
+ // Point p6(rootCS, {0_m, 0_m, 6_m});
+ // Point p7(rootCS, {0_m, 0_m, 7_m});
+ // Point p8(rootCS, {0_m, 0_m, 8_m});
+ // Point p9(rootCS, {0_m, 0_m, 9_m});
+ Point p10(rootCS, {0_m, 0_m, 10_m});
+
+ // get a unit vector
+ Vector<dimensionless_d> v1(rootCS, {0, 0, 1});
+
+ // // get a geometrical path of points
+ // Path P1({p0,p1,p2,p3,p4,p5,p6,p7,p8,p9,p10});
+
+ // construct a Straight Propagator given the uniform refractive index environment
+ StraightPropagator SP(env);
+
+ // store the outcome of the Propagate method to paths_
+ auto const paths_ = SP.propagate(p0, p10, 1_m);
+
+ // perform checks to paths_ components
+ for (auto const& path : paths_) {
+ CHECK((path.time_ / 1_s) - ((34_m / (3 * constants::c)) / 1_s) ==
+ Approx(0).margin(absMargin));
+ CHECK(path.average_refractivity_ == Approx(1));
+ CHECK(path.emit_.getComponents() == v1.getComponents());
+ CHECK(path.receive_.getComponents() == v1.getComponents());
+ // CHECK(std::equal(P1.begin(), P1.end(), path.points_.begin(),[]
+ // (Point a, Point b) { return (a - b).norm() / 1_m < 1e-5;}));
+ //TODO:THINK ABOUT THE POINTS IN THE SIGNALPATH.H
+ }
+
+ CHECK(paths_.size() == 1);
+ }
+
+ SECTION("Straight Propagator w/ Exponential Refractive Index") {
+
+// logging::set_level(logging::level::info);
+// corsika_logger->set_pattern("[%n:%^%-8l%$] custom pattern: %v");
+
+ // create an environment with exponential refractive index (n_0 = 1 & lambda = 0)
+ using ExpoRIndex = ExponentialRefractiveIndex<HomogeneousMedium
+ <IRefractiveIndexModel<IMediumModel>>>;
+
+ using EnvType = Environment<IRefractiveIndexModel<IMediumModel>>;
+ EnvType env1;
+
+ //get another coordinate system
+ const CoordinateSystemPtr rootCS1 = env1.getCoordinateSystem();
+
+ auto Medium1 = EnvType::createNode<Sphere>(
+ Point{rootCS1, 0_m, 0_m, 0_m}, 1_km * std::numeric_limits<double>::infinity());
+
+ auto const props1 =
+ Medium1
+ ->setModelProperties<ExpoRIndex>( 1, 0 / 1_m,
+ 1_kg / (1_m * 1_m * 1_m),
+ NuclearComposition(
+ std::vector<Code>{Code::Nitrogen},
+ std::vector<float>{1.f}));
+
+ env1.getUniverse()->addChild(std::move(Medium1));
+
+ // get some points
+ Point pp0(rootCS1, {0_m, 0_m, 0_m});
+// Point pp1(rootCS1, {0_m, 0_m, 1_m});
+// Point pp2(rootCS1, {0_m, 0_m, 2_m});
+// Point pp3(rootCS1, {0_m, 0_m, 3_m});
+// Point pp4(rootCS1, {0_m, 0_m, 4_m});
+// Point pp5(rootCS1, {0_m, 0_m, 5_m});
+// Point pp6(rootCS1, {0_m, 0_m, 6_m});
+// Point pp7(rootCS1, {0_m, 0_m, 7_m});
+// Point pp8(rootCS1, {0_m, 0_m, 8_m});
+// Point pp9(rootCS1, {0_m, 0_m, 9_m});
+ Point pp10(rootCS1, {0_m, 0_m, 10_m});
+
+ // get a unit vector
+ Vector<dimensionless_d> vv1(rootCS1, {0, 0, 1});
+
+ // construct a Straight Propagator given the exponential refractive index environment
+ StraightPropagator SP1(env1);
+
+ // store the outcome of Propagate method to paths1_
+ auto const paths1_ = SP1.propagate(pp0, pp10, 1_m);
+
+ // perform checks to paths1_ components (this is just a sketch for now)
+ for (auto const& path :paths1_) {
+ CHECK( (path.time_ / 1_s) - ((34_m / (3 * constants::c)) / 1_s)
+ == Approx(0).margin(absMargin) );
+ CHECK( path.average_refractivity_ == Approx(1) );
+ CHECK( path.emit_.getComponents() == vv1.getComponents() );
+ CHECK( path.receive_.getComponents() == vv1.getComponents() );
+ }
+
+ CHECK( paths1_.size() == 1 );
+
+ /*
+ * A second environment with another exponential refractive index
+ */
+
+ // create an environment with exponential refractive index (n_0 = 2 & lambda = 2)
+ using ExpoRIndex = ExponentialRefractiveIndex<HomogeneousMedium
+ <IRefractiveIndexModel<IMediumModel>>>;
+
+ using EnvType = Environment<IRefractiveIndexModel<IMediumModel>>;
+ EnvType env2;
+
+ //get another coordinate system
+ const CoordinateSystemPtr rootCS2 = env2.getCoordinateSystem();
+
+ auto Medium2 = EnvType::createNode<Sphere>(
+ Point{rootCS2, 0_m, 0_m, 0_m}, 1_km * std::numeric_limits<double>::infinity());
+
+ auto const props2 =
+ Medium2
+ ->setModelProperties<ExpoRIndex>( 2, 2 / 1_m,
+ 1_kg / (1_m * 1_m * 1_m),
+ NuclearComposition(
+ std::vector<Code>{Code::Nitrogen},
+ std::vector<float>{1.f}));
+
+ env2.getUniverse()->addChild(std::move(Medium2));
+
+ // get some points
+ Point ppp0(rootCS2, {0_m, 0_m, 0_m});
+ Point ppp10(rootCS2, {0_m, 0_m, 10_m});
+
+ // get a unit vector
+ Vector<dimensionless_d> vvv1(rootCS2, {0, 0, 1});
+
+ // construct a Straight Propagator given the exponential refractive index environment
+ StraightPropagator SP2(env2);
+
+ // store the outcome of Propagate method to paths1_
+ auto const paths2_ = SP2.propagate(ppp0, ppp10, 1_m);
+
+ // perform checks to paths1_ components (this is just a sketch for now)
+ for (auto const& path :paths2_) {
+ CHECK( (path.time_ / 1_s) - ((3.177511688_m / (3 * constants::c)) / 1_s)
+ == Approx(0).margin(absMargin) );
+ CHECK( path.average_refractivity_ == Approx(0.210275935) );
+ CHECK( path.emit_.getComponents() == vvv1.getComponents() );
+ CHECK( path.receive_.getComponents() == vvv1.getComponents() );
+ }
+
+ CHECK( paths2_.size() == 1 );
+
+ }
+
+// SECTION("Construct a ZHS process.") {
+//
+// // TODO: construct the environment for the propagator
+//
+// // here we just use a deque to store the antenna
+// std::deque<TimeDomainAntenna> antennas;
+//
+// // TODO: add time domain antennas to the deque
+//
+// // and now create ZHS with this antenna collection
+// // and with a straight propagator
+// ZHS<decltype(antennas), StraightPropagator> zhs(antennas, env);
+// // TODO: do we need the explicit type declaration?
+//
+// // here is where the simulation happens
+//
+// // and call zhs.writeOutput() at the end.
+// }
+
+} // END: TEST_CASE("Radio", "[processes]")