From a6b865b5f1ebb9be09b9ed2dfad88c3bcab90e15 Mon Sep 17 00:00:00 2001 From: Nikos Karastathis <n.karastathis@kit.edu> Date: Wed, 6 Jul 2022 17:41:03 +0200 Subject: [PATCH] add getters for sampling rate and start time --- corsika/detail/modules/radio/CoREAS.inl | 10 ++-- corsika/detail/modules/radio/RadioProcess.inl | 2 +- corsika/detail/modules/radio/ZHS.inl | 46 +++++++++---------- .../radio/antennas/TimeDomainAntenna.inl | 4 ++ .../radio/antennas/TimeDomainAntenna.hpp | 21 +++++++-- 5 files changed, 49 insertions(+), 34 deletions(-) diff --git a/corsika/detail/modules/radio/CoREAS.inl b/corsika/detail/modules/radio/CoREAS.inl index 1322ef163..cf0793c84 100644 --- a/corsika/detail/modules/radio/CoREAS.inl +++ b/corsika/detail/modules/radio/CoREAS.inl @@ -195,7 +195,7 @@ namespace corsika { // CoREAS calculation -> get ElectricFieldVector for "midPoint" ElectricFieldVector EVmid_ = (path.emit_.cross(path.emit_.cross(beta_))) / midDoppler_ / path.R_distance_ * constants_ * - antenna.sample_rate_; + antenna.getSampleRate(); ElectricFieldVector EV1_{EVmid_}; ElectricFieldVector EV2_{EVmid_ * (-1.0)}; @@ -214,7 +214,7 @@ namespace corsika { endPointReceiveTime_ = midPointReceiveTime_ - 0.5 * deltaT_; } - TimeType const gridResolution_{1 / antenna.sample_rate_}; + TimeType const gridResolution_{1 / antenna.getSampleRate()}; deltaT_ = endPointReceiveTime_ - startPointReceiveTime_; // redistribute contributions over time scale defined by the observation @@ -308,19 +308,19 @@ namespace corsika { // calculate electric field vector for startpoint ElectricFieldVector EV1_ = (paths1[i].emit_.cross(paths1[i].emit_.cross(beta_))) / preDoppler_ / - paths1[i].R_distance_ * constants_ * antenna.sample_rate_; + paths1[i].R_distance_ * constants_ * antenna.getSampleRate(); // calculate electric field vector for endpoint ElectricFieldVector EV2_ = (paths2[i].emit_.cross(paths2[i].emit_.cross(beta_))) / postDoppler_ / - paths2[i].R_distance_ * constants_ * (-1.0) * antenna.sample_rate_; + paths2[i].R_distance_ * constants_ * (-1.0) * antenna.getSampleRate(); if ((preDoppler_ < 1.e-9) || (postDoppler_ < 1.e-9)) { // CORSIKA_LOG_ERROR("Doppler factors are less than 1.e-9 // for this track"); - TimeType const gridResolution_{1 / antenna.sample_rate_}; + TimeType const gridResolution_{1 / antenna.getSampleRate()}; TimeType deltaT_{endPointReceiveTime_ - startPointReceiveTime_}; if (abs(deltaT_) < (gridResolution_)) { diff --git a/corsika/detail/modules/radio/RadioProcess.inl b/corsika/detail/modules/radio/RadioProcess.inl index 8198a3665..f53aba72d 100644 --- a/corsika/detail/modules/radio/RadioProcess.inl +++ b/corsika/detail/modules/radio/RadioProcess.inl @@ -83,7 +83,7 @@ namespace corsika { // before the next event for (auto& antenna : antennas_.getAntennas()) { antenna.endOfShower(event_, this->implementation().algorithm, - antenna.sample_rate_ * 1_s); + antenna.getSampleRate() * 1_s); antenna.reset(); } diff --git a/corsika/detail/modules/radio/ZHS.inl b/corsika/detail/modules/radio/ZHS.inl index 29f175e0a..eb381c09f 100644 --- a/corsika/detail/modules/radio/ZHS.inl +++ b/corsika/detail/modules/radio/ZHS.inl @@ -47,7 +47,7 @@ namespace corsika { for (size_t i{0}; i < midPaths.size(); i++) { double const uTimesK{beta.dot(midPaths[i].emit_) / betaModule}; double const sinTheta2{1. - uTimesK * uTimesK}; - LengthType const lambda{constants::c / antenna.sample_rate_}; + LengthType const lambda{constants::c / antenna.getSampleRate()}; double const fraunhLimit{sinTheta2 * trackLength * trackLength / midPaths[i].R_distance_ / lambda * 2 * M_PI}; // Checks if we are in fraunhoffer domain @@ -92,9 +92,9 @@ namespace corsika { } // end if statement for time structure double const startBin{std::floor( - (detectionTime1 - antenna.start_time_) * antenna.sample_rate_ + 0.5)}; + (detectionTime1 - antenna.getStartTime()) * antenna.getSampleRate() + 0.5)}; double const endBin{std::floor( - (detectionTime2 - antenna.start_time_) * antenna.sample_rate_ + 0.5)}; + (detectionTime2 - antenna.getStartTime()) * antenna.getSampleRate() + 0.5)}; auto const betaPerp{ newMidPaths[k].emit_.cross(beta.cross(newMidPaths[k].emit_))}; @@ -104,14 +104,14 @@ namespace corsika { // track contained in bin // if not in Cerenkov angle then if (std::fabs(denominator) > 1.e-15) { - double const f{std::fabs((detectionTime2 * antenna.sample_rate_ - - detectionTime1 * antenna.sample_rate_))}; + double const f{std::fabs((detectionTime2 * antenna.getSampleRate() - + detectionTime1 * antenna.getSampleRate()))}; VectorPotential const Vp = betaPerp * sign * constants * f / denominator / newMidPaths[k].R_distance_; antenna.receive(detectionTime2, betaPerp, Vp); } else { // If emission in Cerenkov angle => approximation - double const f{time2 * antenna.sample_rate_ - - time1 * antenna.sample_rate_}; + double const f{time2 * antenna.getSampleRate() - + time1 * antenna.getSampleRate()}; VectorPotential const Vp = betaPerp * sign * constants * f / newMidPaths[k].R_distance_; antenna.receive(detectionTime2, betaPerp, Vp); @@ -121,8 +121,8 @@ namespace corsika { int const numberOfBins{static_cast<int>(endBin - startBin)}; // first contribution/ plus 1 bin minus 0.5 from new antenna ruonding double f{std::fabs(startBin + 0.5 - - (detectionTime1 - antenna.start_time_) * - antenna.sample_rate_)}; + (detectionTime1 - antenna.getStartTime()) * + antenna.getSampleRate())}; VectorPotential Vp = betaPerp * sign * constants * f / denominator / newMidPaths[k].R_distance_; antenna.receive(detectionTime1, betaPerp, Vp); @@ -130,12 +130,12 @@ namespace corsika { for (int it{1}; it < numberOfBins; ++it) { Vp = betaPerp * constants / denominator / newMidPaths[k].R_distance_; antenna.receive( - detectionTime1 + static_cast<double>(it) / antenna.sample_rate_, + detectionTime1 + static_cast<double>(it) / antenna.getSampleRate(), betaPerp, Vp); } // end loop over bins in which potential vector is not zero // final contribution// f +0.5 from new antenna rounding - f = std::fabs((detectionTime2 - antenna.start_time_) * - antenna.sample_rate_ + + f = std::fabs((detectionTime2 - antenna.getStartTime()) * + antenna.getSampleRate() + 0.5 - endBin); Vp = betaPerp * sign * constants * f / denominator / newMidPaths[k].R_distance_; @@ -171,9 +171,9 @@ namespace corsika { } // end if statement for time structure double const startBin{std::floor( - (detectionTime1 - antenna.start_time_) * antenna.sample_rate_ + 0.5)}; + (detectionTime1 - antenna.getStartTime()) * antenna.getSampleRate() + 0.5)}; double const endBin{std::floor( - (detectionTime2 - antenna.start_time_) * antenna.sample_rate_ + 0.5)}; + (detectionTime2 - antenna.getStartTime()) * antenna.getSampleRate() + 0.5)}; auto const betaPerp{midPaths[i].emit_.cross(beta.cross(midPaths[i].emit_))}; double const denominator{1. - @@ -183,15 +183,15 @@ namespace corsika { // track contained in bin // if not in Cerenkov angle then if (std::fabs(denominator) > 1.e-15) { - double const f{std::fabs((detectionTime2 * antenna.sample_rate_ - - detectionTime1 * antenna.sample_rate_))}; + double const f{std::fabs((detectionTime2 * antenna.getSampleRate() - + detectionTime1 * antenna.getSampleRate()))}; VectorPotential const Vp = betaPerp * sign * constants * f / denominator / midPaths[i].R_distance_; antenna.receive(detectionTime2, betaPerp, Vp); } else { // If emission in Cerenkov angle => approximation - double const f{endTime * antenna.sample_rate_ - - startTime * antenna.sample_rate_}; + double const f{endTime * antenna.getSampleRate() - + startTime * antenna.getSampleRate()}; VectorPotential const Vp = betaPerp * sign * constants * f / midPaths[i].R_distance_; antenna.receive(detectionTime2, betaPerp, Vp); @@ -202,8 +202,8 @@ namespace corsika { // TODO: should we check for Cerenkov angle? // first contribution double f{std::fabs(startBin + 0.5 - - (detectionTime1 - antenna.start_time_) * - antenna.sample_rate_)}; + (detectionTime1 - antenna.getStartTime()) * + antenna.getSampleRate())}; VectorPotential Vp = betaPerp * sign * constants * f / denominator / midPaths[i].R_distance_; antenna.receive(detectionTime1, betaPerp, Vp); @@ -211,12 +211,12 @@ namespace corsika { for (int it{1}; it < numberOfBins; ++it) { Vp = betaPerp * sign * constants / denominator / midPaths[i].R_distance_; antenna.receive( - detectionTime1 + static_cast<double>(it) / antenna.sample_rate_, + detectionTime1 + static_cast<double>(it) / antenna.getSampleRate(), betaPerp, Vp); } // end loop over bins in which potential vector is not zero // final contribution - f = std::fabs((detectionTime2 - antenna.start_time_) * - antenna.sample_rate_ + + f = std::fabs((detectionTime2 - antenna.getStartTime()) * + antenna.getSampleRate() + 0.5 - endBin); Vp = betaPerp * sign * constants * f / denominator / midPaths[i].R_distance_; diff --git a/corsika/detail/modules/radio/antennas/TimeDomainAntenna.inl b/corsika/detail/modules/radio/antennas/TimeDomainAntenna.inl index e8448eeea..b0a7e4b55 100644 --- a/corsika/detail/modules/radio/antennas/TimeDomainAntenna.inl +++ b/corsika/detail/modules/radio/antennas/TimeDomainAntenna.inl @@ -99,6 +99,10 @@ namespace corsika { inline auto const& TimeDomainAntenna::getAxis() const { return time_axis_; } + inline InverseTimeType const& TimeDomainAntenna::getSampleRate() const { return sample_rate_; } + + inline TimeType const& TimeDomainAntenna::getStartTime() const { return start_time_; } + inline void TimeDomainAntenna::reset() { std::fill(waveformEX_.begin(), waveformEX_.end(), 0); std::fill(waveformEY_.begin(), waveformEY_.end(), 0); diff --git a/corsika/modules/radio/antennas/TimeDomainAntenna.hpp b/corsika/modules/radio/antennas/TimeDomainAntenna.hpp index ee411c29d..88bf1dc46 100644 --- a/corsika/modules/radio/antennas/TimeDomainAntenna.hpp +++ b/corsika/modules/radio/antennas/TimeDomainAntenna.hpp @@ -21,11 +21,6 @@ namespace corsika { */ class TimeDomainAntenna : public Antenna<TimeDomainAntenna> { - public: - // import the methods from the antenna - - // label this as a time-domain antenna. - static constexpr bool is_time_domain{true}; TimeType const start_time_; ///< The start time of this waveform. TimeType const duration_; ///< The duration of this waveform. @@ -37,6 +32,12 @@ namespace corsika { TimeType const ground_hit_time_; ///< The time the primary particle hits the ground. std::vector<long double> const time_axis_; ///< The time axis corresponding to the electric field. + public: + // import the methods from the antenna + + // label this as a time-domain antenna. + static constexpr bool is_time_domain{true}; + using Antenna<TimeDomainAntenna>::getName; using Antenna<TimeDomainAntenna>::getLocation; @@ -111,6 +112,16 @@ namespace corsika { */ auto const& getAxis() const; + /** + * Returns the sampling rate of the time domain antenna. + */ + InverseTimeType const& getSampleRate() const; + + /** + * Returns the start time of detection for the time domain antenna. + */ + TimeType const& getStartTime() const; + /** * Reset the antenna before starting a new simulation. */ -- GitLab