From 8bf837b85ae68f3bf8b9e1cda0964b34550041e3 Mon Sep 17 00:00:00 2001 From: ralfulrich <ralf.ulrich@kit.edu> Date: Fri, 5 Feb 2021 10:57:38 +0100 Subject: [PATCH] Max's comments, and a bug... --- corsika/detail/modules/ObservationPlane.inl | 8 +- corsika/detail/modules/tracking/Intersect.inl | 16 +-- .../tracking/TrackingLeapFrogCurved.inl | 125 +++++++++--------- .../tracking/TrackingLeapFrogStraight.inl | 32 ++--- .../modules/tracking/TrackingStraight.inl | 18 +-- .../tracking/TrackingLeapFrogCurved.hpp | 15 ++- corsika/modules/tracking/TrackingStraight.hpp | 15 +-- do-clang-format.py | 1 + src/modules/qgsjetII/code_generator.py | 4 +- 9 files changed, 108 insertions(+), 126 deletions(-) diff --git a/corsika/detail/modules/ObservationPlane.inl b/corsika/detail/modules/ObservationPlane.inl index f3e77a4c8..33cb556ec 100644 --- a/corsika/detail/modules/ObservationPlane.inl +++ b/corsika/detail/modules/ObservationPlane.inl @@ -59,13 +59,9 @@ namespace corsika { corsika::setup::Stack::particle_type const& particle, corsika::setup::Trajectory const& trajectory) { - auto const& volumeNode = particle.getNode(); - typedef typename std::remove_const_t< - std::remove_reference_t<decltype(volumeNode->getModelProperties())>> - medium_type; Intersections const intersection = - setup::Tracking::intersect<corsika::setup::Stack::particle_type, medium_type>( - particle, plane_, volumeNode->getModelProperties()); + setup::Tracking::intersect<corsika::setup::Stack::particle_type>(particle, + plane_); TimeType const timeOfIntersection = intersection.getEntry(); CORSIKA_LOG_TRACE("particle={}, pos={}, dir={}, plane={}, timeOfIntersection={}", particle.asString(), particle.getPosition(), diff --git a/corsika/detail/modules/tracking/Intersect.inl b/corsika/detail/modules/tracking/Intersect.inl index fef16a955..79cfa6984 100644 --- a/corsika/detail/modules/tracking/Intersect.inl +++ b/corsika/detail/modules/tracking/Intersect.inl @@ -38,7 +38,7 @@ namespace corsika { // first check, where we leave the current volume // this assumes our convention, that all Volume primitives must be convex // thus, the last entry is always the exit point - const Intersections time_intersections_curr = + Intersections const time_intersections_curr = TDerived::intersect(particle, volumeNode); CORSIKA_LOG_TRACE("curr node {}, parent node {}, hasIntersections={} ", fmt::ptr(&volumeNode), fmt::ptr(volumeNode.getParent()), @@ -56,16 +56,16 @@ namespace corsika { // where do we collide with any of the next-tree-level volumes // entirely contained by currentLogicalVolumeNode - for (const auto& node : volumeNode.getChildNodes()) { + for (auto const& node : volumeNode.getChildNodes()) { - const Intersections time_intersections = TDerived::intersect(particle, *node); + Intersections const time_intersections = TDerived::intersect(particle, *node); if (!time_intersections.hasIntersections()) { continue; } CORSIKA_LOG_DEBUG("intersection times with child volume {} : enter {} s, exit {} s", fmt::ptr(node), time_intersections.getEntry() / 1_s, time_intersections.getExit() / 1_s); - const auto t_entry = time_intersections.getEntry(); - const auto t_exit = time_intersections.getExit(); + auto const t_entry = time_intersections.getEntry(); + auto const t_exit = time_intersections.getExit(); CORSIKA_LOG_TRACE("children t-entry: {}, t-exit: {}, smaller? {} ", t_entry, t_exit, t_entry <= minTime); // note, theoretically t can even be smaller than 0 since we @@ -86,14 +86,14 @@ namespace corsika { // current volume for (node_type* node : volumeNode.getExcludedNodes()) { - const Intersections time_intersections = TDerived::intersect(particle, *node); + Intersections const time_intersections = TDerived::intersect(particle, *node); if (!time_intersections.hasIntersections()) { continue; } CORSIKA_LOG_DEBUG( "intersection times with exclusion volume {} : enter {} s, exit {} s", fmt::ptr(node), time_intersections.getEntry() / 1_s, time_intersections.getExit() / 1_s); - const auto t_entry = time_intersections.getEntry(); - const auto t_exit = time_intersections.getExit(); + auto const t_entry = time_intersections.getEntry(); + auto const t_exit = time_intersections.getExit(); CORSIKA_LOG_TRACE("children t-entry: {}, t-exit: {}, smaller? {} ", t_entry, t_exit, t_entry <= minTime); // note, theoretically t can even be smaller than 0 since we diff --git a/corsika/detail/modules/tracking/TrackingLeapFrogCurved.inl b/corsika/detail/modules/tracking/TrackingLeapFrogCurved.inl index 18a1bf0d2..8dc9f1f3a 100644 --- a/corsika/detail/modules/tracking/TrackingLeapFrogCurved.inl +++ b/corsika/detail/modules/tracking/TrackingLeapFrogCurved.inl @@ -103,9 +103,9 @@ namespace corsika { LengthType const gyroradius = (pAlongB_delta * 1_V / (constants::c * abs(chargeNumber) * magnitudeB * 1_eV)); - const double maxRadians = 0.01; - const LengthType steplimit = 2 * cos(maxRadians) * sin(maxRadians) * gyroradius; - const TimeType steplimit_time = steplimit / initialVelocity.getNorm(); + double const maxRadians = 0.01; + LengthType const steplimit = 2 * cos(maxRadians) * sin(maxRadians) * gyroradius; + TimeType const steplimit_time = steplimit / initialVelocity.getNorm(); CORSIKA_LOG_DEBUG("gyroradius {}, steplimit: {} = {}", gyroradius, steplimit, steplimit_time); @@ -113,7 +113,7 @@ namespace corsika { // intersection auto [minTime, minNode] = nextIntersect(particle, steplimit_time); - const auto k = + auto const k = chargeNumber * constants::cSquared * 1_eV / (particle.getEnergy() * 1_V); return std::make_tuple( LeapFrogTrajectory(position, initialVelocity, magneticfield, k, @@ -121,10 +121,9 @@ namespace corsika { minNode); // next volume node } - template <typename TParticle, typename TMedium> - inline Intersections Tracking::intersect(const TParticle& particle, - const Sphere& sphere, - const TMedium& medium) { + template <typename TParticle> + inline Intersections Tracking::intersect(TParticle const& particle, + Sphere const& sphere) { if (sphere.getRadius() == 1_km * std::numeric_limits<double>::infinity()) { return Intersections(); @@ -132,7 +131,9 @@ namespace corsika { int const chargeNumber = particle.getChargeNumber(); auto const& position = particle.getPosition(); - MagneticFieldVector const& magneticfield = medium.getMagneticField(position); + auto const* currentLogicalVolumeNode = particle.getNode(); + MagneticFieldVector const& magneticfield = + currentLogicalVolumeNode->getModelProperties().getMagneticField(position); VelocityVector const velocity = particle.getMomentum() / particle.getEnergy() * constants::c; @@ -144,26 +145,23 @@ namespace corsika { bool const isParallel = (projectedDirectionSqrNorm == 0 * square(1_T)); if (chargeNumber == 0 || magneticfield.getNorm() == 0_T || isParallel) { - return tracking_line::Tracking::intersect<TParticle, TMedium>(particle, sphere, - medium); + return tracking_line::Tracking::intersect<TParticle>(particle, sphere); } bool const numericallyInside = sphere.contains(particle.getPosition()); - const auto absVelocity = velocity.getNorm(); + auto const absVelocity = velocity.getNorm(); auto energy = particle.getEnergy(); auto k = chargeNumber * constants::cSquared * 1_eV / (absVelocity * energy * 1_V); - auto const denom = (directionBefore.cross(magneticfield)).getSquaredNorm() * k * k; - const double a = - ((directionBefore.cross(magneticfield)).dot(position - sphere.getCenter()) * k + - 1) * - 4 / (1_m * 1_m * denom); - const double b = directionBefore.dot(position - sphere.getCenter()) * 8 / - (denom * 1_m * 1_m * 1_m); - const double c = ((position - sphere.getCenter()).getSquaredNorm() - - (sphere.getRadius() * sphere.getRadius())) * - 4 / (denom * 1_m * 1_m * 1_m * 1_m); + auto const direction_B_perp = directionBefore.cross(magneticfield); + auto const denom = direction_B_perp.getSquaredNorm() * k * k; + Vector<length_d> const deltaPos = position - sphere.getCenter(); + double const a = (direction_B_perp.dot(deltaPos) * k + 1) * 4 / (1_m * 1_m * denom); + double const b = directionBefore.dot(deltaPos) * 8 / (denom * 1_m * 1_m * 1_m); + double const c = + (deltaPos.getSquaredNorm() - (sphere.getRadius() * sphere.getRadius())) * 4 / + (denom * 1_m * 1_m * 1_m * 1_m); CORSIKA_LOG_TRACE("denom={}, a={}, b={}, c={}", denom, a, b, c); std::complex<double>* solutions = quartic_solver::solve_quartic(0, a, b, c); LengthType d_enter, d_exit; @@ -225,9 +223,9 @@ namespace corsika { return Intersections(d_enter / absVelocity, d_exit / absVelocity); } - template <typename TParticle, typename TMedium> + template <typename TParticle> inline Intersections Tracking::intersect(TParticle const& particle, - Plane const& plane, TMedium const& medium) { + Plane const& plane) { int chargeNumber; if (is_nucleus(particle.getPID())) { @@ -236,9 +234,6 @@ namespace corsika { chargeNumber = get_charge_number(particle.getPID()); } auto const* currentLogicalVolumeNode = particle.getNode(); - auto magneticfield = - currentLogicalVolumeNode->getModelProperties().getMagneticField( - particle.getPosition()); VelocityVector const velocity = particle.getMomentum() / particle.getEnergy() * constants::c; auto const absVelocity = velocity.getNorm(); @@ -246,70 +241,68 @@ namespace corsika { velocity.normalized(); // determine steplength to next volume Point const position = particle.getPosition(); - if (chargeNumber != 0 && - abs(plane.getNormal().dot(velocity.cross(magneticfield))) * 1_s / 1_m / 1_T > - 1e-6) { + auto const magneticfield = + currentLogicalVolumeNode->getModelProperties().getMagneticField(position); + + if (chargeNumber != 0 && abs(plane.getNormal().dot(velocity.cross(magneticfield))) > + 1e-6_T * 1_m / 1_s) { auto const* currentLogicalVolumeNode = particle.getNode(); - auto magneticfield = - currentLogicalVolumeNode->getModelProperties().getMagneticField( - particle.getPosition()); - auto k = - chargeNumber * constants::c * 1_eV / (particle.getMomentum().getNorm() * 1_V); - - if (direction.dot(plane.getNormal()) * direction.dot(plane.getNormal()) - - (plane.getNormal().dot(position - plane.getCenter()) * - plane.getNormal().dot(direction.cross(magneticfield)) * 2 * k) < - 0) { + auto const magneticfield = + currentLogicalVolumeNode->getModelProperties().getMagneticField(position); + auto const k = + chargeNumber * (constants::c * 1_eV / 1_V) / particle.getMomentum().getNorm(); + + auto const direction_B_perp = direction.cross(magneticfield); + auto const denom = plane.getNormal().dot(direction_B_perp) * k; + auto const sqrtArg = + direction.dot(plane.getNormal()) * direction.dot(plane.getNormal()) - + (plane.getNormal().dot(position - plane.getCenter()) * denom * 2); + + if (sqrtArg < 0) { return Intersections(std::numeric_limits<double>::infinity() * 1_s); } + double const sqrSqrtArg = sqrt(sqrtArg); + auto const norm_projected = + direction.dot(plane.getNormal()) / direction.getNorm(); + LengthType const MaxStepLength1 = (sqrSqrtArg - norm_projected) / denom; + LengthType const MaxStepLength2 = (-sqrSqrtArg - norm_projected) / denom; - LengthType const MaxStepLength1 = - (sqrt(direction.dot(plane.getNormal()) * direction.dot(plane.getNormal()) - - (plane.getNormal().dot(position - plane.getCenter()) * - plane.getNormal().dot(direction.cross(magneticfield)) * 2 * k)) - - direction.dot(plane.getNormal()) / direction.getNorm()) / - (plane.getNormal().dot(direction.cross(magneticfield)) * k); - - LengthType const MaxStepLength2 = - (-sqrt(direction.dot(plane.getNormal()) * direction.dot(plane.getNormal()) - - (plane.getNormal().dot(position - plane.getCenter()) * - plane.getNormal().dot(direction.cross(magneticfield)) * 2 * k)) - - direction.dot(plane.getNormal()) / direction.getNorm()) / - (plane.getNormal().dot(direction.cross(magneticfield)) * k); - + // check: both intersections in past if (MaxStepLength1 <= 0_m && MaxStepLength2 <= 0_m) { return Intersections(std::numeric_limits<double>::infinity() * 1_s); + + // check: next intersection is MaxStepLength2 } else if (MaxStepLength1 <= 0_m || MaxStepLength2 < MaxStepLength1) { CORSIKA_LOG_TRACE(" steplength to obs plane 2: {} ", MaxStepLength2); return Intersections( MaxStepLength2 * - (direction + direction.cross(magneticfield) * MaxStepLength2 * k / 2) - .getNorm() / + (direction + direction_B_perp * MaxStepLength2 * k / 2).getNorm() / absVelocity); + + // check: next intersections is MaxStepLength1 } else if (MaxStepLength2 <= 0_m || MaxStepLength1 < MaxStepLength2) { CORSIKA_LOG_TRACE(" steplength to obs plane 2: {} ", MaxStepLength1); return Intersections( MaxStepLength1 * - (direction + direction.cross(magneticfield) * MaxStepLength2 * k / 2) - .getNorm() / + (direction + direction_B_perp * MaxStepLength1 * k / 2).getNorm() / absVelocity); } - CORSIKA_LOG_WARN("Particle wasn't tracked with curved trajectory -> straight"); + CORSIKA_LOG_WARN( + "Particle wasn't tracked with curved trajectory -> straight (is this an " + "ERROR?)"); } // end if curved-tracking - return tracking_line::Tracking::intersect(particle, plane, medium); + return tracking_line::Tracking::intersect(particle, plane); } template <typename TParticle, typename TBaseNodeType> - inline Intersections Tracking::intersect(const TParticle& particle, - const TBaseNodeType& volumeNode) { - const Sphere* sphere = dynamic_cast<const Sphere*>(&volumeNode.getVolume()); - if (sphere) { - return intersect(particle, *sphere, volumeNode.getModelProperties()); - } + inline Intersections Tracking::intersect(TParticle const& particle, + TBaseNodeType const& volumeNode) { + Sphere const* sphere = dynamic_cast<Sphere const*>(&volumeNode.getVolume()); + if (sphere) { return intersect(particle, *sphere); } throw std::runtime_error( "The Volume type provided is not supported in intersect(particle, node)"); } diff --git a/corsika/detail/modules/tracking/TrackingLeapFrogStraight.inl b/corsika/detail/modules/tracking/TrackingLeapFrogStraight.inl index 7f913bd17..61ad3793d 100644 --- a/corsika/detail/modules/tracking/TrackingLeapFrogStraight.inl +++ b/corsika/detail/modules/tracking/TrackingLeapFrogStraight.inl @@ -28,7 +28,7 @@ namespace corsika { VelocityVector initialVelocity = particle.getMomentum() / particle.getEnergy() * constants::c; - const Point initialPosition = particle.getPosition(); + Point const initialPosition = particle.getPosition(); CORSIKA_LOG_DEBUG( "TrackingB pid: {}" " , E = {} GeV", @@ -53,17 +53,17 @@ namespace corsika { const int chargeNumber = particle.getChargeNumber(); auto magneticfield = volumeNode->getModelProperties().getMagneticField(initialPosition); - const auto magnitudeB = magneticfield.getNorm(); + auto const magnitudeB = magneticfield.getNorm(); CORSIKA_LOG_DEBUG("field={} uT, chargeNumber={}, magnitudeB={} uT", magneticfield.getComponents() / 1_uT, chargeNumber, magnitudeB / 1_T); bool const no_deflection = chargeNumber == 0 || magnitudeB == 0_T; // check, where the first halve-step direction has geometric intersections - const auto [initialTrack, initialTrackNextVolume] = + auto const [initialTrack, initialTrackNextVolume] = tracking_line::Tracking::getTrack(particle); { [[maybe_unused]] auto& initialTrackNextVolume_dum = initialTrackNextVolume; } - const auto initialTrackLength = initialTrack.getLength(1); + auto const initialTrackLength = initialTrack.getLength(1); CORSIKA_LOG_DEBUG("initialTrack(0)={}, initialTrack(1)={}, initialTrackLength={}", initialTrack.getPosition(0).getCoordinates(), @@ -94,8 +94,8 @@ namespace corsika { // need to follow strongly curved trajectories segment-wise, // at least if we don't employ concepts as "Helix // Trajectories" or similar - const double maxRadians = 0.01; - const LengthType steplimit = 2 * cos(maxRadians) * sin(maxRadians) * gyroradius; + double const maxRadians = 0.01; + LengthType const steplimit = 2 * cos(maxRadians) * sin(maxRadians) * gyroradius; CORSIKA_LOG_DEBUG("gyroradius {}, Steplimit: {}", gyroradius, steplimit); // calculate first halve step for "steplimit" @@ -110,12 +110,12 @@ namespace corsika { CORSIKA_LOG_DEBUG("first halve step length {}, steplimit={}, initialTrackLength={}", firstHalveSteplength, steplimit, initialTrackLength); // perform the first halve-step - const Point position_mid = initialPosition + direction * firstHalveSteplength; - const auto k = - chargeNumber * constants::c * 1_eV / (particle.getMomentum().getNorm() * 1_V); - const auto new_direction = + Point const position_mid = initialPosition + direction * firstHalveSteplength; + auto const k = + chargeNumber * (constants::c * 1_eV / 1_V) / particle.getMomentum().getNorm(); + auto const new_direction = direction + direction.cross(magneticfield) * firstHalveSteplength * 2 * k; - const auto new_direction_norm = new_direction.getNorm(); // by design this is >1 + auto const new_direction_norm = new_direction.getNorm(); // by design this is >1 CORSIKA_LOG_DEBUG( "position_mid={}, new_direction={}, (new_direction_norm)={}, deflection={}", position_mid.getCoordinates(), new_direction.getComponents(), @@ -126,7 +126,7 @@ namespace corsika { // check, where the second halve-step direction has geometric intersections particle.setPosition(position_mid); particle.setMomentum(new_direction * absMomentum); - const auto [finalTrack, finalTrackNextVolume] = + auto const [finalTrack, finalTrackNextVolume] = tracking_line::Tracking::getTrack(particle); particle.setPosition(initialPosition); // this is not nice... particle.setMomentum(initialMomentum); // this is not nice... @@ -160,12 +160,12 @@ namespace corsika { // perform the second halve-step auto const new_direction_normalized = new_direction.normalized(); - const Point finalPosition = + Point const finalPosition = position_mid + new_direction_normalized * secondHalveStepLength; - const LengthType totalStep = firstHalveSteplength + secondHalveStepLength; - const auto delta_pos = finalPosition - initialPosition; - const auto distance = delta_pos.getNorm(); + LengthType const totalStep = firstHalveSteplength + secondHalveStepLength; + auto const delta_pos = finalPosition - initialPosition; + auto const distance = delta_pos.getNorm(); return std::make_tuple( StraightTrajectory( diff --git a/corsika/detail/modules/tracking/TrackingStraight.inl b/corsika/detail/modules/tracking/TrackingStraight.inl index 7fc57b21b..411ab54ff 100644 --- a/corsika/detail/modules/tracking/TrackingStraight.inl +++ b/corsika/detail/modules/tracking/TrackingStraight.inl @@ -48,9 +48,9 @@ namespace corsika::tracking_line { minNode); // next volume node } - template <typename TParticle, typename TMedium> + template <typename TParticle> inline Intersections Tracking::intersect(TParticle const& particle, - Sphere const& sphere, TMedium const&) { + Sphere const& sphere) { auto const delta = particle.getPosition() - sphere.getCenter(); auto const velocity = particle.getMomentum() / particle.getEnergy() * constants::c; auto const vSqNorm = velocity.getSquaredNorm(); @@ -73,20 +73,14 @@ namespace corsika::tracking_line { inline Intersections Tracking::intersect(TParticle const& particle, TBaseNodeType const& volumeNode) { Sphere const* sphere = dynamic_cast<Sphere const*>(&volumeNode.getVolume()); - if (sphere) { - typedef typename std::remove_const_t< - std::remove_reference_t<decltype(volumeNode.getModelProperties())>> - medium_type; - return Tracking::intersect<TParticle, medium_type>(particle, *sphere, - volumeNode.getModelProperties()); - } + if (sphere) { return Tracking::intersect<TParticle>(particle, *sphere); } throw std::runtime_error( "The Volume type provided is not supported in Intersect(particle, node)"); } - template <typename TParticle, typename TMedium> - inline Intersections Tracking::intersect(TParticle const& particle, Plane const& plane, - TMedium const&) { + template <typename TParticle> + inline Intersections Tracking::intersect(TParticle const& particle, + Plane const& plane) { auto const delta = plane.getCenter() - particle.getPosition(); auto const velocity = particle.getMomentum() / particle.getEnergy() * constants::c; auto const n = plane.getNormal(); diff --git a/corsika/modules/tracking/TrackingLeapFrogCurved.hpp b/corsika/modules/tracking/TrackingLeapFrogCurved.hpp index 05073691e..5eb281608 100644 --- a/corsika/modules/tracking/TrackingLeapFrogCurved.hpp +++ b/corsika/modules/tracking/TrackingLeapFrogCurved.hpp @@ -60,17 +60,18 @@ namespace corsika { template <typename TParticle> auto getTrack(TParticle const& particle); - template <typename TParticle, typename TMedium> - static Intersections intersect(TParticle const& particle, Sphere const& sphere, - TMedium const& medium); + //! find intersection of Sphere with Track + template <typename TParticle> + static Intersections intersect(TParticle const& particle, Sphere const& sphere); + //! find intersection of Volume node with Track of particle template <typename TParticle, typename TBaseNodeType> static Intersections intersect(TParticle const& particle, - TBaseNodeType const& volumeNode); + TBaseNodeType const& node); - template <typename TParticle, typename TMedium> - static Intersections intersect(TParticle const& particle, Plane const& plane, - TMedium const& medium); + //! find intersection of Plane with Track + template <typename TParticle> + static Intersections intersect(TParticle const& particle, Plane const& plane); protected: /** diff --git a/corsika/modules/tracking/TrackingStraight.hpp b/corsika/modules/tracking/TrackingStraight.hpp index b77f8fd21..4b498bfa5 100644 --- a/corsika/modules/tracking/TrackingStraight.hpp +++ b/corsika/modules/tracking/TrackingStraight.hpp @@ -41,19 +41,16 @@ namespace corsika::tracking_line { auto getTrack(TParticle const& particle); //! find intersection of Sphere with Track - template <typename TParticle, typename TMedium> - static Intersections intersect(TParticle const& particle, Sphere const& sphere, - TMedium const&); + template <typename TParticle> + static Intersections intersect(TParticle const& particle, Sphere const& sphere); - //! find intersection of Volume with Track + //! find intersection of Volume node with Track of particle template <typename TParticle, typename TBaseNodeType> - static Intersections intersect(TParticle const& particle, - TBaseNodeType const& volumeNode); + static Intersections intersect(TParticle const& particle, TBaseNodeType const& node); //! find intersection of Plane with Track - template <typename TParticle, typename TMedium> - static Intersections intersect(TParticle const& particle, Plane const& plane, - TMedium const&); + template <typename TParticle> + static Intersections intersect(TParticle const& particle, Plane const& plane); }; } // namespace corsika::tracking_line diff --git a/do-clang-format.py b/do-clang-format.py index 24973f669..3bff0e619 100755 --- a/do-clang-format.py +++ b/do-clang-format.py @@ -83,6 +83,7 @@ cmd += " -style=file" version = subp.check_output(cmd.split() + ["--version"]).decode("utf-8") print (version) +print ("Note: the clang-format version has an impact on the result. Make sure you are consistent with current CI. Consider \'--docker\' option.") if args.apply: for filename in filelist: diff --git a/src/modules/qgsjetII/code_generator.py b/src/modules/qgsjetII/code_generator.py index ce7fcae26..8321aad19 100755 --- a/src/modules/qgsjetII/code_generator.py +++ b/src/modules/qgsjetII/code_generator.py @@ -93,7 +93,7 @@ def read_qgsjetII_codes(filename, particle_db): if len(line)==0 or line[0] == '#': continue line = line.split('#')[0] - print (line) + print ('QGSJetII codes: ', line) identifier, model_code, xsType = line.split() try: particle_db[identifier]["qgsjetII_code"] = int(model_code) @@ -202,7 +202,7 @@ if __name__ == "__main__": particle_db = load_particledb(sys.argv[1]) read_qgsjetII_codes(sys.argv[2], particle_db) set_default_qgsjetII_definition(particle_db) - + with open("Generated.inc", "w") as f: print("// this file is automatically generated\n// edit at your own risk!\n", file=f) print(generate_qgsjetII_start(), file=f) -- GitLab