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 runAndEvent/RE04/src/RE04Trajectory.cc
0027 /// \brief Implementation of the RE04Trajectory class
0028 //
0029 //
0030 #include "RE04Trajectory.hh"
0031
0032 #include "RE04TrajectoryPoint.hh"
0033
0034 #include "G4AttDef.hh"
0035 #include "G4AttDefStore.hh"
0036 #include "G4AttValue.hh"
0037 #include "G4ParticleTable.hh"
0038 #include "G4UIcommand.hh"
0039 #include "G4UnitsTable.hh"
0040
0041 // #define G4ATTDEBUG
0042 #ifdef G4ATTDEBUG
0043 #  include "G4AttCheck.hh"
0044 #endif
0045
0046 G4ThreadLocal G4Allocator<RE04Trajectory>* faTrajAllocator = 0;
0047
0048 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0049 RE04Trajectory::RE04Trajectory()
0050   : G4VTrajectory(),
0051     fPositionRecord(0),
0052     fTrackID(0),
0053     fParentID(0),
0054     fPDGEncoding(0),
0055     fPDGCharge(0.0),
0056     fParticleName(""),
0057     fInitialKineticEnergy(0.),
0058     fInitialMomentum(G4ThreeVector())
0059 {
0060   ;
0061 }
0062
0063 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0064 RE04Trajectory::RE04Trajectory(const G4Track* aTrack)
0065 {
0066   G4ParticleDefinition* fpParticleDefinition = aTrack->GetDefinition();
0067   fParticleName = fpParticleDefinition->GetParticleName();
0068   fPDGCharge = fpParticleDefinition->GetPDGCharge();
0069   fPDGEncoding = fpParticleDefinition->GetPDGEncoding();
0070   fTrackID = aTrack->GetTrackID();
0071   fParentID = aTrack->GetParentID();
0072   fInitialKineticEnergy = aTrack->GetKineticEnergy();
0073   fInitialMomentum = aTrack->GetMomentum();
0074   fPositionRecord = new TrajectoryPointContainer();
0075   // Following is for the first trajectory point
0076   fPositionRecord->push_back(new RE04TrajectoryPoint(aTrack->GetPosition(), aTrack->GetMaterial()));
0077 }
0078
0079 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0080 RE04Trajectory::RE04Trajectory(RE04Trajectory& right) : G4VTrajectory()
0081 {
0082   fParticleName = right.fParticleName;
0083   fPDGCharge = right.fPDGCharge;
0084   fPDGEncoding = right.fPDGEncoding;
0085   fTrackID = right.fTrackID;
0086   fParentID = right.fParentID;
0087   fInitialKineticEnergy = right.fInitialKineticEnergy;
0088   fInitialMomentum = right.fInitialMomentum;
0089   fPositionRecord = new TrajectoryPointContainer();
0090
0091   for (size_t i = 0; i < right.fPositionRecord->size(); i++) {
0092     RE04TrajectoryPoint* rightPoint = (RE04TrajectoryPoint*)((*(right.fPositionRecord))[i]);
0093     fPositionRecord->push_back(new RE04TrajectoryPoint(*rightPoint));
0094   }
0095 }
0096
0097 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0098 RE04Trajectory::~RE04Trajectory()
0099 {
0100   if (fPositionRecord) {
0101     //  fPositionRecord->clearAndDestroy();
0102     size_t i;
0103     for (i = 0; i < fPositionRecord->size(); i++) {
0104       delete (*fPositionRecord)[i];
0105     }
0106     fPositionRecord->clear();
0107     delete fPositionRecord;
0108   }
0109 }
0110
0111 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0112 void RE04Trajectory::ShowTrajectory(std::ostream& os) const
0113 {
0114   // Invoke the default implementation in G4VTrajectory...
0115   G4VTrajectory::ShowTrajectory(os);
0116   // ... or override with your own code here.
0117 }
0118
0119 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0120 void RE04Trajectory::DrawTrajectory() const
0121 {
0122   // Invoke the default implementation in G4VTrajectory...
0123   G4VTrajectory::DrawTrajectory();
0124   // ... or override with your own code here.
0125 }
0126
0127 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0128 const std::map<G4String, G4AttDef>* RE04Trajectory::GetAttDefs() const
0129 {
0130   G4bool isNew;
0131   std::map<G4String, G4AttDef>* store = G4AttDefStore::GetInstance("RE04Trajectory", isNew);
0132   if (isNew) {
0133     G4String id("ID");
0134     (*store)[id] = G4AttDef(id, "Track ID", "Physics", "", "G4int");
0135
0136     G4String pid("PID");
0137     (*store)[pid] = G4AttDef(pid, "Parent ID", "Physics", "", "G4int");
0138
0139     G4String pn("PN");
0140     (*store)[pn] = G4AttDef(pn, "Particle Name", "Physics", "", "G4String");
0141
0142     G4String ch("Ch");
0143     (*store)[ch] = G4AttDef(ch, "Charge", "Physics", "e+", "G4double");
0144
0145     G4String pdg("PDG");
0146     (*store)[pdg] = G4AttDef(pdg, "PDG Encoding", "Physics", "", "G4int");
0147
0148     G4String ike("IKE");
0149     (*store)[ike] = G4AttDef(ike, "Initial kinetic energy", "Physics", "G4BestUnit", "G4double");
0150
0151     G4String iMom("IMom");
0152     (*store)[iMom] = G4AttDef(iMom, "Initial momentum", "Physics", "G4BestUnit", "G4ThreeVector");
0153
0154     G4String iMag("IMag");
0155     (*store)[iMag] =
0156       G4AttDef(iMag, "Initial momentum magnitude", "Physics", "G4BestUnit", "G4double");
0157
0158     G4String ntp("NTP");
0159     (*store)[ntp] = G4AttDef(ntp, "No. of points", "Physics", "", "G4int");
0160   }
0161   return store;
0162 }
0163
0164 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0165 std::vector<G4AttValue>* RE04Trajectory::CreateAttValues() const
0166 {
0167   std::vector<G4AttValue>* values = new std::vector<G4AttValue>;
0168
0169   values->push_back(G4AttValue("ID", G4UIcommand::ConvertToString(fTrackID), ""));
0170
0171   values->push_back(G4AttValue("PID", G4UIcommand::ConvertToString(fParentID), ""));
0172
0173   values->push_back(G4AttValue("PN", fParticleName, ""));
0174
0175   values->push_back(G4AttValue("Ch", G4UIcommand::ConvertToString(fPDGCharge), ""));
0176
0177   values->push_back(G4AttValue("PDG", G4UIcommand::ConvertToString(fPDGEncoding), ""));
0178
0179   values->push_back(G4AttValue("IKE", G4BestUnit(fInitialKineticEnergy, "Energy"), ""));
0180
0181   values->push_back(G4AttValue("IMom", G4BestUnit(fInitialMomentum, "Energy"), ""));
0182
0183   values->push_back(G4AttValue("IMag", G4BestUnit(fInitialMomentum.mag(), "Energy"), ""));
0184
0185   values->push_back(G4AttValue("NTP", G4UIcommand::ConvertToString(GetPointEntries()), ""));
0186
0187 #ifdef G4ATTDEBUG
0188   G4cout << G4AttCheck(values, GetAttDefs());
0189 #endif
0190
0191   return values;
0192 }
0193
0194 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0195 void RE04Trajectory::AppendStep(const G4Step* aStep)
0196 {
0197   fPositionRecord->push_back(new RE04TrajectoryPoint(aStep->GetPostStepPoint()->GetPosition(),
0198                                                      aStep->GetPreStepPoint()->GetMaterial()));
0199 }
0200
0201 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0202 G4ParticleDefinition* RE04Trajectory::GetParticleDefinition()
0203 {
0204   return (G4ParticleTable::GetParticleTable()->FindParticle(fParticleName));
0205 }
0206
0207 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0208 void RE04Trajectory::MergeTrajectory(G4VTrajectory* secondTrajectory)
0209 {
0210   if (!secondTrajectory) return;
0211
0212   RE04Trajectory* seco = (RE04Trajectory*)secondTrajectory;
0213   G4int ent = seco->GetPointEntries();
0214   for (G4int i = 1; i < ent; i++)  // initial point of the second trajectory
0215                                    // should not be merged
0216   {
0217     fPositionRecord->push_back((*(seco->fPositionRecord))[i]);
0218     //    fPositionRecord->push_back(seco->fPositionRecord->removeAt(1));
0219   }
0220   delete (*seco->fPositionRecord)[0];
0221   seco->fPositionRecord->clear();
0222 }