RE04/src/RE04Trajectory.cc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 /// \file RE04Trajectory.cc
0027 /// \brief Implementation of the RE04Trajectory class
0028
0029 #include "RE04Trajectory.hh"
0030
0031 #include "RE04TrajectoryPoint.hh"
0032
0033 #include "G4AttDef.hh"
0034 #include "G4AttDefStore.hh"
0035 #include "G4AttValue.hh"
0036 #include "G4ParticleTable.hh"
0037 #include "G4UIcommand.hh"
0038 #include "G4UnitsTable.hh"
0039
0040 // #define G4ATTDEBUG
0041 #ifdef G4ATTDEBUG
0042 #  include "G4AttCheck.hh"
0043 #endif
0044
0045 G4ThreadLocal G4Allocator<RE04Trajectory>* faTrajAllocator = 0;
0046
0047 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0048 RE04Trajectory::RE04Trajectory()
0049   : G4VTrajectory(),
0050     fPositionRecord(0),
0051     fTrackID(0),
0052     fParentID(0),
0053     fPDGEncoding(0),
0054     fPDGCharge(0.0),
0055     fParticleName(""),
0056     fInitialKineticEnergy(0.),
0057     fInitialMomentum(G4ThreeVector())
0058 {
0059   ;
0060 }
0061
0062 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0063 RE04Trajectory::RE04Trajectory(const G4Track* aTrack)
0064 {
0065   G4ParticleDefinition* fpParticleDefinition = aTrack->GetDefinition();
0066   fParticleName = fpParticleDefinition->GetParticleName();
0067   fPDGCharge = fpParticleDefinition->GetPDGCharge();
0068   fPDGEncoding = fpParticleDefinition->GetPDGEncoding();
0069   fTrackID = aTrack->GetTrackID();
0070   fParentID = aTrack->GetParentID();
0071   fInitialKineticEnergy = aTrack->GetKineticEnergy();
0072   fInitialMomentum = aTrack->GetMomentum();
0073   fPositionRecord = new TrajectoryPointContainer();
0074   // Following is for the first trajectory point
0075   fPositionRecord->push_back(new RE04TrajectoryPoint(aTrack->GetPosition(), aTrack->GetMaterial()));
0076 }
0077
0078 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0079 RE04Trajectory::RE04Trajectory(RE04Trajectory& right) : G4VTrajectory()
0080 {
0081   fParticleName = right.fParticleName;
0082   fPDGCharge = right.fPDGCharge;
0083   fPDGEncoding = right.fPDGEncoding;
0084   fTrackID = right.fTrackID;
0085   fParentID = right.fParentID;
0086   fInitialKineticEnergy = right.fInitialKineticEnergy;
0087   fInitialMomentum = right.fInitialMomentum;
0088   fPositionRecord = new TrajectoryPointContainer();
0089
0090   for (size_t i = 0; i < right.fPositionRecord->size(); i++) {
0091     RE04TrajectoryPoint* rightPoint = (RE04TrajectoryPoint*)((*(right.fPositionRecord))[i]);
0092     fPositionRecord->push_back(new RE04TrajectoryPoint(*rightPoint));
0093   }
0094 }
0095
0096 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0097 RE04Trajectory::~RE04Trajectory()
0098 {
0099   if (fPositionRecord) {
0100     //  fPositionRecord->clearAndDestroy();
0101     size_t i;
0102     for (i = 0; i < fPositionRecord->size(); i++) {
0103       delete (*fPositionRecord)[i];
0104     }
0105     fPositionRecord->clear();
0106     delete fPositionRecord;
0107   }
0108 }
0109
0110 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0111 void RE04Trajectory::ShowTrajectory(std::ostream& os) const
0112 {
0113   // Invoke the default implementation in G4VTrajectory...
0114   G4VTrajectory::ShowTrajectory(os);
0115   // ... or override with your own code here.
0116 }
0117
0118 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0119 void RE04Trajectory::DrawTrajectory() const
0120 {
0121   // Invoke the default implementation in G4VTrajectory...
0122   G4VTrajectory::DrawTrajectory();
0123   // ... or override with your own code here.
0124 }
0125
0126 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0127 const std::map<G4String, G4AttDef>* RE04Trajectory::GetAttDefs() const
0128 {
0129   G4bool isNew;
0130   std::map<G4String, G4AttDef>* store = G4AttDefStore::GetInstance("RE04Trajectory", isNew);
0131   if (isNew) {
0132     G4String id("ID");
0133     (*store)[id] = G4AttDef(id, "Track ID", "Physics", "", "G4int");
0134
0135     G4String pid("PID");
0136     (*store)[pid] = G4AttDef(pid, "Parent ID", "Physics", "", "G4int");
0137
0138     G4String pn("PN");
0139     (*store)[pn] = G4AttDef(pn, "Particle Name", "Physics", "", "G4String");
0140
0141     G4String ch("Ch");
0142     (*store)[ch] = G4AttDef(ch, "Charge", "Physics", "e+", "G4double");
0143
0144     G4String pdg("PDG");
0145     (*store)[pdg] = G4AttDef(pdg, "PDG Encoding", "Physics", "", "G4int");
0146
0147     G4String ike("IKE");
0148     (*store)[ike] = G4AttDef(ike, "Initial kinetic energy", "Physics", "G4BestUnit", "G4double");
0149
0150     G4String iMom("IMom");
0151     (*store)[iMom] = G4AttDef(iMom, "Initial momentum", "Physics", "G4BestUnit", "G4ThreeVector");
0152
0153     G4String iMag("IMag");
0154     (*store)[iMag] =
0155       G4AttDef(iMag, "Initial momentum magnitude", "Physics", "G4BestUnit", "G4double");
0156
0157     G4String ntp("NTP");
0158     (*store)[ntp] = G4AttDef(ntp, "No. of points", "Physics", "", "G4int");
0159   }
0160   return store;
0161 }
0162
0163 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0164 std::vector<G4AttValue>* RE04Trajectory::CreateAttValues() const
0165 {
0166   std::vector<G4AttValue>* values = new std::vector<G4AttValue>;
0167
0168   values->push_back(G4AttValue("ID", G4UIcommand::ConvertToString(fTrackID), ""));
0169
0170   values->push_back(G4AttValue("PID", G4UIcommand::ConvertToString(fParentID), ""));
0171
0172   values->push_back(G4AttValue("PN", fParticleName, ""));
0173
0174   values->push_back(G4AttValue("Ch", G4UIcommand::ConvertToString(fPDGCharge), ""));
0175
0176   values->push_back(G4AttValue("PDG", G4UIcommand::ConvertToString(fPDGEncoding), ""));
0177
0178   values->push_back(G4AttValue("IKE", G4BestUnit(fInitialKineticEnergy, "Energy"), ""));
0179
0180   values->push_back(G4AttValue("IMom", G4BestUnit(fInitialMomentum, "Energy"), ""));
0181
0182   values->push_back(G4AttValue("IMag", G4BestUnit(fInitialMomentum.mag(), "Energy"), ""));
0183
0184   values->push_back(G4AttValue("NTP", G4UIcommand::ConvertToString(GetPointEntries()), ""));
0185
0186 #ifdef G4ATTDEBUG
0187   G4cout << G4AttCheck(values, GetAttDefs());
0188 #endif
0189
0190   return values;
0191 }
0192
0193 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0194 void RE04Trajectory::AppendStep(const G4Step* aStep)
0195 {
0196   fPositionRecord->push_back(new RE04TrajectoryPoint(aStep->GetPostStepPoint()->GetPosition(),
0197                                                      aStep->GetPreStepPoint()->GetMaterial()));
0198 }
0199
0200 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0201 G4ParticleDefinition* RE04Trajectory::GetParticleDefinition()
0202 {
0203   return (G4ParticleTable::GetParticleTable()->FindParticle(fParticleName));
0204 }
0205
0206 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0207 void RE04Trajectory::MergeTrajectory(G4VTrajectory* secondTrajectory)
0208 {
0209   if (!secondTrajectory) return;
0210
0211   RE04Trajectory* seco = (RE04Trajectory*)secondTrajectory;
0212   G4int ent = seco->GetPointEntries();
0213   for (G4int i = 1; i < ent; i++)  // initial point of the second trajectory
0214                                    // should not be merged
0215   {
0216     fPositionRecord->push_back((*(seco->fPositionRecord))[i]);
0217     //    fPositionRecord->push_back(seco->fPositionRecord->removeAt(1));
0218   }
0219   delete (*seco->fPositionRecord)[0];
0220   seco->fPositionRecord->clear();
0221 }