fixed getMaxStepLength in radio process - radio showers compile now with some tweaking

corsika/modules/radio/CoREAS.hpp

@@ -346,30 +346,6 @@ namespace corsika {
       } // End of looping over antennas
     } // End of simulate method
 
-    /**
-     * 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 {
-
-      // TODO : This is where we control the maximum step size
-      // of a particle track in order to maintain the accuracy
-      // of the particular formalism.
-      //
-      // This is part of the ZHS / CoReas formalisms and can
-      // be related from the magnetic field / acceleration, charge,
-      // etc. of the particle.
-      return 1000000000000_m;
-    }
-
-  }; // END: class RadioProcess
+  }; // END: class CoREAS
 
 } // namespace corsika
\ No newline at end of file

corsika/modules/radio/RadioProcess.hpp

@@ -133,6 +133,21 @@ namespace corsika {
 
       }
 
+    /**
+   * 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.
+   */
+    LengthType getMaxStepLength(setup::Stack::particle_type const& vParticle,
+                                setup::Trajectory const& vTrack) const {
+      return meter * std::numeric_limits<double>::infinity();
+    }
+
   }; // END: class RadioProcess
 
 } // namespace corsika

examples/radio_shower.cpp

@@ -99,12 +99,15 @@ void registerRandomStreams(const int seed) {
 
 template <typename TInterface>
 using MyExtraEnv =
-    UniformRefractiveIndex<MediumPropertyModel<UniformMagneticField<TInterface>>>;
+UniformRefractiveIndex<MediumPropertyModel<UniformMagneticField<TInterface>>>;
+
+//template <typename T>
+//using MyExtraEnv = MediumPropertyModel<UniformMagneticField<T>>;
 
 int main(int argc, char** argv) {
 
   corsika_logger->set_pattern("[%n:%^%-8l%$] %s:%#: %v");
-  logging::set_level(logging::level::trace);
+  logging::set_level(logging::level::info);
 
   CORSIKA_LOG_INFO("vertical_EAS");
 
@@ -120,10 +123,14 @@ int main(int argc, char** argv) {
   // initialize random number sequence(s)
   registerRandomStreams(seed);
 
-  // setup environment
-  using EnvironmentInterface =
-  IRefractiveIndexModel<IMediumPropertyModel<IMagneticFieldModel<IMediumModel>>>;
-  using EnvType = Environment<EnvironmentInterface>;
+  // setup 2 environments (use only one)
+
+//  using EnvironmentInterface =
+//  IRefractiveIndexModel<IMediumPropertyModel<IMagneticFieldModel<IMediumModel>>>;
+//  using EnvType = Environment<EnvironmentInterface>;
+//  EnvType env;
+
+  using EnvType = setup::Environment;
   EnvType env;
   CoordinateSystemPtr const& rootCS = env.getCoordinateSystem();
   Point const center{rootCS, 0_m, 0_m, 0_m};
@@ -131,18 +138,25 @@ int main(int argc, char** argv) {
   const auto point1{Point(env.getCoordinateSystem(), 50_m, 50_m, 50_m)};
   const auto point2{Point(env.getCoordinateSystem(), 25_m, 25_m, 25_m)};
   // the antennas
-  TimeDomainAntenna ant1("antenna1", point1, 0_s, 100_s, 1/1e-6_s);
-  TimeDomainAntenna ant2("antenna2", point2, 0_s, 100_s, 1/1e-6_s);
+  TimeDomainAntenna ant1("antenna1", point1, 0_s, 1_s, 1/1e-6_s);
+  TimeDomainAntenna ant2("antenna2", point2, 0_s, 1_s, 1/1e-6_s);
   // the detector
   AntennaCollection<TimeDomainAntenna> detector;
   detector.addAntenna(ant1);
   detector.addAntenna(ant2);
   auto builder = make_layered_spherical_atmosphere_builder<
-      EnvironmentInterface, MyExtraEnv>::create(center,
+      setup::EnvironmentInterface, MyExtraEnv>::create(center,
                                                        constants::EarthRadius::Mean, 1.000327,
                                                        Medium::AirDry1Atm,
                                                        MagneticFieldVector{rootCS, 0_T,
                                                                            50_uT, 0_T});
+  // builder with refractive index interface
+//  auto builder = make_layered_spherical_atmosphere_builder<
+//      EnvironmentInterface, MyExtraEnv>::create(center,
+//                                                       constants::EarthRadius::Mean, 1.000327,
+//                                                       Medium::AirDry1Atm,
+//                                                       MagneticFieldVector{rootCS, 0_T,
+//                                                                           50_uT, 0_T});
 
   builder.setNuclearComposition(
       {{Code::Nitrogen, Code::Oxygen},

tests/modules/testRadio.cpp

@@ -164,7 +164,7 @@ TEST_CASE("Radio", "[processes]") {
     // create times for the antenna
     const TimeType t1{0_s}; // TODO: initialization of times to antennas! particle hits the observation level should be zero
     const TimeType t2{100_s};
-    const InverseTimeType t3{1/1_s};
+    const InverseTimeType t3{1/1e+2_s};
     const TimeType t4{11_s};
 
     // check that I can create an antenna at (1, 2, 3)
@@ -172,6 +172,7 @@ TEST_CASE("Radio", "[processes]") {
     TimeDomainAntenna ant2("antenna_name2", point2, t1, t2, t3);
 //    TimeDomainAntenna ant3("antenna1", point1, 0_s, 2_s, 1/1e-7_s);
 
+//    std::cout << "static cast " << static_cast<int>(1/1000) << std::endl;
 
     // construct a radio detector instance to store our antennas
     AntennaCollection<TimeDomainAntenna> detector;