min distance is no longer steplength

Framework/Cascade/Cascade.h

@@ -408,7 +408,7 @@ namespace corsika::cascade {
 
       // determine steplength for the magnetic field
       // because Steplength should not be min_distance
-      /*
+      
       int chargeNumber;
 	    if (corsika::particles::IsNucleus(vParticle.GetPID())) {
 	      chargeNumber = vParticle.GetNuclearZ();
@@ -419,16 +419,21 @@ namespace corsika::cascade {
 	    auto magneticfield = currentLogicalVolumeNode->GetModelProperties().GetMagneticField(vParticle.GetPosition());
       geometry::Vector<dimensionless_d> const directionBefore = vParticle.GetMomentum().normalized();
 	    auto c = directionBefore.cross(magneticfield) * chargeNumber * corsika::units::constants::c * 1_eV / 
-               (vParticle.GetMomentum().norm() * 1_V)
+               (vParticle.GetMomentum().norm() * 1_V);
 	    LengthType Steplength = min_distance;
 	    if (chargeNumber != 0) {
 	      Steplength = sqrt(2 / c.squaredNorm() * (sqrt(c.squaredNorm() * min_distance * min_distance + 1) -1));
         std::cout << "Steplength " << Steplength << std::endl;
-	    } 
-	    */
+	    }
+      if (Steplength == 0_m) {
+        Steplength = min_distance;
+      }
+	    
 	    // This formula hasnt been tested
 	    
-      auto [position, direction, L2] = fTracking.MagneticStep(vParticle, min_distance);
+      //auto [position, direction, L2] = fTracking.MagneticStep(vParticle, min_distance);
+      auto [position, direction, L2] = fTracking.MagneticStep(vParticle, Steplength);
+      
       //histL2(L2);
       histLlog2(L2);
       int pdg = static_cast<int>(particles::GetPDG(vParticle.GetPID()));
@@ -481,7 +486,7 @@ namespace corsika::cascade {
               
               
               
-     int chargeNumber = 0;
+     //int chargeNumber = 0;
         if (corsika::particles::IsNucleus(vParticle.GetPID())) {
         	chargeNumber = vParticle.GetNuclearZ();
         } else {

Processes/ObservationPlane/ObservationPlane.cc

@@ -71,7 +71,7 @@ LengthType ObservationPlane::MaxStepLength(setup::Stack::ParticleType const& vPa
   auto magneticfield = currentLogicalVolumeNode->GetModelProperties().GetMagneticField(vParticle.GetPosition());
   auto direction = trajectory.GetV0().normalized();
   
-  if (chargeNumber != 0 && plane_.GetNormal().dot(trajectory.GetV0().cross(magneticfield)) * 1_s / 1_m / 1_T > 1e-6) {
+  if (chargeNumber != 0 && abs(plane_.GetNormal().dot(trajectory.GetV0().cross(magneticfield))) * 1_s / 1_m / 1_T > 1e-6) {
     auto const* currentLogicalVolumeNode = vParticle.GetNode();
     auto magneticfield = currentLogicalVolumeNode->GetModelProperties().GetMagneticField(vParticle.GetPosition());
     auto k = chargeNumber * corsika::units::constants::c * 1_eV / (vParticle.GetMomentum().norm() * 1_V);