← Back to branch summary

~chris-rogers/maus/emr_mc_digitization

« back to all changes in this revision

Viewing changes to src/common_cpp/Optics/LinearApproximationOpticsModel.cc

  • Committer: Chris Rogers
  • Date: 2014-04-16 11:48:45 UTC
  • mfrom: (707 merge)
  • mto: This revision was merged to the branch mainline in revision 711.
  • Revision ID: chris.rogers@stfc.ac.uk-20140416114845-h3u3q7pdcxkxvovs
Update to trunk

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/common_cpp/Optics/LinearApproximationOpticsModel.cc
23
23
#include <vector>
24
24
#include <iostream>
25
25
 
 
26
#include "CLHEP/Units/PhysicalConstants.h"
 
27
 
26
28
#include "src/common_cpp/Optics/CovarianceMatrix.hh"
27
29
#include "src/common_cpp/Optics/PhaseSpaceVector.hh"
28
30
#include "src/common_cpp/Recon/Global/Particle.hh"
29
 
#include "src/common_cpp/Recon/Global/TrackPoint.hh"
30
31
#include "Maths/Vector.hh"
31
32
 
32
 
#include "Interface/Squeal.hh"
 
33
#include "Utils/Exception.hh"
33
34
 
34
35
namespace MAUS {
35
36
 
36
37
const TransferMap * LinearApproximationOpticsModel::CalculateTransferMap(
37
 
    const std::vector<recon::global::TrackPoint> & start_plane_hits,
38
 
    const std::vector<recon::global::TrackPoint> & station_hits)
 
38
    const std::vector<MAUS::PhaseSpaceVector> & start_plane_hits,
 
39
    const std::vector<MAUS::PhaseSpaceVector> & station_hits)
39
40
    const {
40
41
  const MAUS::DataStructure::Global::PID particle_id
41
 
    = MAUS::DataStructure::Global::PID(start_plane_hits[0].particle_id());
 
42
    = MAUS::DataStructure::Global::PID(reference_primary_.GetParticleId());
 
43
  const double mass
 
44
    = recon::global::Particle::GetInstance().GetMass(particle_id);
42
45
 
43
 
  double start_plane = start_plane_hits[0].z();
 
46
  double start_plane = reference_primary_.GetPosition().z();
44
47
 
45
48
  double hit_total = 0.0;
46
 
  std::vector<recon::global::TrackPoint>::const_iterator hit;
 
49
  std::vector<MAUS::PhaseSpaceVector>::const_iterator hit;
47
50
  for (hit = station_hits.begin(); hit != station_hits.end(); ++hit) {
48
 
    hit_total += hit->z();
 
51
    const double energy = hit->energy();
 
52
    const double momentum = ::sqrt(energy*energy - mass*mass);
 
53
    const double beta = momentum / energy;
 
54
    const double delta_t = hit->time() - time_offset_
 
55
                         - reference_primary_.GetTime();
 
56
    const double delta_z = beta * ::CLHEP::c_light * delta_t;
 
57
    hit_total += (reference_primary_.GetPosition().z() + delta_z);
49
58
  }
50
59
  double end_plane = hit_total / station_hits.size();
51
60
 
52
 
  const double mass
53
 
    = recon::global::Particle::GetInstance()->GetMass(particle_id);
54
 
 
55
61
  return new LinearApproximationTransferMap(start_plane, end_plane, mass);
56
62
}
57
63
 
78
84
 
79
85
  double energy = vector.E();
80
86
  if (mass_ > energy) {
81
 
    throw(Squeal(Squeal::nonRecoverable,
 
87
    throw(Exception(Exception::nonRecoverable,
82
88
                 "Attempting to transport a particle that is off mass shell.",
83
89
                 "MAUS::LinearApproximationOpticsModel::Transport()"));
84
90
  }