Line data    Source code

       1             : #include "crpropa/ParticleState.h"
       2             : #include "crpropa/Units.h"
       3             : #include "crpropa/Common.h"
       4             : #include "crpropa/ParticleID.h"
       5             : #include "crpropa/ParticleMass.h"
       6             : 
       7             : #include "HepPID/ParticleIDMethods.hh"
       8             : 
       9             : #include <cstdlib>
      10             : #include <sstream>
      11             : 
      12             : namespace crpropa {
      13             : 
      14    16859556 : ParticleState::ParticleState(int id, double E, Vector3d pos, Vector3d dir): id(0), energy(0.), position(0.), direction(0.), pmass(0.), charge(0.)
      15             : {
      16    16859556 :         setId(id);
      17    16859556 :         setEnergy(E);
      18    16859556 :         setPosition(pos);
      19    16859556 :         setDirection(dir);
      20    16859556 : }
      21             : 
      22    24161173 : void ParticleState::setPosition(const Vector3d &pos) {
      23             :         position = pos;
      24    24161173 : }
      25             : 
      26    12548182 : const Vector3d &ParticleState::getPosition() const {
      27    12548182 :         return position;
      28             : }
      29             : 
      30    17344794 : void ParticleState::setDirection(const Vector3d &dir) {
      31             :         direction = dir / dir.getR();
      32    17344794 : }
      33             : 
      34      525404 : const Vector3d &ParticleState::getDirection() const {
      35      525404 :         return direction;
      36             : }
      37             : 
      38    20366015 : void ParticleState::setEnergy(double newEnergy) {
      39    20366015 :         energy = std::max(0., newEnergy); // prevent negative energies
      40    20366015 : }
      41             : 
      42      653619 : double ParticleState::getEnergy() const {
      43      653619 :         return energy;
      44             : }
      45             : 
      46           1 : double ParticleState::getRigidity() const {
      47           1 :         return fabs(energy / charge);
      48             : }
      49             : 
      50    20332793 : void ParticleState::setId(int newId) {
      51    20332793 :         id = newId;
      52    20332793 :         if (isNucleus(id)) {
      53      146401 :                 pmass = nuclearMass(id);
      54      146401 :                 charge = chargeNumber(id) * eplus;
      55      146401 :                 if (id < 0)
      56           4 :                         charge *= -1; // anti-nucleus
      57             :         } else {
      58    20186392 :                 if (abs(id) == 11)
      59       16250 :                         pmass = mass_electron;
      60    20186392 :                 charge = HepPID::charge(id) * eplus;
      61             :         }
      62    20332793 : }
      63             : 
      64      565472 : int ParticleState::getId() const {
      65      565472 :         return id;
      66             : }
      67             : 
      68           4 : double ParticleState::getMass() const {
      69           4 :         return pmass;
      70             : }
      71             : 
      72      960335 : double ParticleState::getCharge() const {
      73      960335 :         return charge;
      74             : }
      75             : 
      76      132796 : double ParticleState::getLorentzFactor() const {
      77      132796 :         return energy / (pmass * c_squared);
      78             : }
      79             : 
      80       65429 : void ParticleState::setLorentzFactor(double lf) {
      81       65429 :         lf = std::max(0., lf); // prevent negative Lorentz factors
      82       65429 :         energy = lf * pmass * c_squared;
      83       65429 : }
      84             : 
      85           1 : Vector3d ParticleState::getVelocity() const {
      86           1 :         return direction * c_light;
      87             : }
      88             : 
      89           2 : Vector3d ParticleState::getMomentum() const {
      90           2 :         return direction * (energy / c_light);
      91             : }
      92             : 
      93           0 : std::string ParticleState::getDescription() const {
      94           0 :         std::stringstream ss;
      95           0 :         ss << "Particle " << id << ", ";
      96           0 :         ss << "E = " << energy / EeV << " EeV, ";
      97           0 :         ss << "x = " << position / Mpc << " Mpc, ";
      98           0 :         ss << "p = " << direction;
      99           0 :         return ss.str();
     100           0 : }
     101             : 
     102             : } // namespace crpropa