TACTICScorer.cc 3.66 KB
Newer Older
1 2
//#define USE_Garfield //only use if compiled with Garfield

Warren's avatar
Warren committed
3 4
#include "TACTICScorer.hh"
#include "G4UnitsTable.hh"
Warren's avatar
Warren committed
5
#include "G4RunManager.hh"
6 7 8 9 10

#ifdef USE_Garfield
#include "GARFDRIFT.h"
#endif

Warren's avatar
Warren committed
11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27
using namespace TACTICScorer;

Gas_Scorer::Gas_Scorer(G4String name,G4int Level,G4double ScorerLength,G4int NumberOfSegments, G4int depth) //what do level and depth do?       
:G4VPrimitiveScorer(name, depth),HCID(-1){
  m_ScorerLength = ScorerLength;
  m_NumberOfSegments = NumberOfSegments;
  m_SegmentLength = m_ScorerLength / m_NumberOfSegments;
  m_Level = Level;
  m_Position = G4ThreeVector(-1000,-1000,-1000);
  m_SegmentNumber = -1;
  m_Index = -1;
}

Gas_Scorer::~Gas_Scorer(){}

G4bool Gas_Scorer::ProcessHits(G4Step* aStep, G4TouchableHistory*){

Warren's avatar
Warren committed
28
  G4double* Infos = new G4double[13];
Warren's avatar
Warren committed
29
  m_Position  = aStep->GetPreStepPoint()->GetPosition();
30

Warren's avatar
Warren committed
31
  Infos[0] = G4RunManager::GetRunManager()->GetCurrentEvent()->GetEventID();
Warren's avatar
Warren committed
32

Warren's avatar
Warren committed
33 34 35
  Infos[1] = aStep->GetTrack()->GetTrackID();

  Infos[2] = aStep->GetTrack()->GetParticleDefinition()->GetAtomicNumber();;
36
  
Warren's avatar
Warren committed
37 38 39
  Infos[3] = aStep->GetPreStepPoint()->GetGlobalTime();
  Infos[4] = aStep->GetPreStepPoint()->GetKineticEnergy();
  Infos[5] = aStep->GetTotalEnergyDeposit() - aStep->GetNonIonizingEnergyDeposit();
Warren's avatar
Warren committed
40 41

  m_SegmentNumber = (int)((m_Position.z() + m_ScorerLength / 2.) / m_SegmentLength ) + 1; //Pad number 
Warren's avatar
Warren committed
42 43 44 45 46
  Infos[6] = m_SegmentNumber;
  Infos[7] = m_Position.z();
  Infos[8] = pow(pow(m_Position.x(),2) + pow(m_Position.y(),2),0.5); //R
  Infos[9] = aStep->GetTrack()->GetVertexPosition()[2];
  Infos[10] = aStep->GetTrack()->GetVertexKineticEnergy();
Warren's avatar
Warren committed
47
  G4ThreeVector p_vec = aStep->GetTrack()->GetVertexMomentumDirection();
Warren's avatar
Warren committed
48 49 50
  Infos[11] = acos(p_vec[2]/pow(pow(p_vec[0],2)+pow(p_vec[1],2)+pow(p_vec[2],2),0.5))/deg; //angle relative to z axis (theta);   
  Infos[12] = aStep->GetTrack()->GetTrackLength();
  
Warren's avatar
Warren committed
51 52
  m_DetectorNumber = aStep->GetPreStepPoint()->GetTouchableHandle()->GetCopyNumber(m_Level);
  m_Index = m_DetectorNumber * 1e3 + m_SegmentNumber * 1e6;
53
  
Warren's avatar
Warren committed
54
  if(isnan(Infos[10])) {
55 56 57 58
    aStep->GetTrack()->SetTrackStatus(fStopAndKill);
    return 0;
  }
            
59 60 61
#ifdef USE_Garfield
  G4ThreeVector delta_Position = aStep->GetDeltaPosition();

Warren's avatar
Warren committed
62
  GARFDRIFT(Infos[5]/eV, Infos[3], m_Position/cm, delta_Position/cm, Infos[8]/cm, Infos[6], Infos[2], m_ScorerLength/cm, m_SegmentLength/cm, Infos[0], Infos[2], Infos[11]);
63 64
#endif
  
Warren's avatar
Warren committed
65 66
  map<G4int, G4double**>::iterator it;
  it= EvtMap->GetMap()->find(m_Index);
67
  if(it!=EvtMap->GetMap()->end()){
Warren's avatar
Warren committed
68
    G4double* dummy = *(it->second);
69
    if(Infos[1]==dummy[1]) Infos[4]+=dummy[4]; //accumulate ionisation energy deposit to get total accross pad
Warren's avatar
Warren committed
70 71 72 73
    delete dummy;
  }
  
  EvtMap->set(m_Index, Infos);
74

Warren's avatar
Warren committed
75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107
  return TRUE;

}

void Gas_Scorer::Initialize(G4HCofThisEvent* HCE){
    EvtMap = new NPS::HitsMap<G4double*>(GetMultiFunctionalDetector()->GetName(), GetName());
    if (HCID < 0) {
        HCID = GetCollectionID(0);
    }
    HCE->AddHitsCollection(HCID, (G4VHitsCollection*)EvtMap);
}

void Gas_Scorer::EndOfEvent(G4HCofThisEvent*){}

void Gas_Scorer::clear(){
    std::map<G4int, G4double**>::iterator    MapIterator;
    for (MapIterator = EvtMap->GetMap()->begin() ; MapIterator != EvtMap->GetMap()->end() ; MapIterator++){
        delete *(MapIterator->second);
    }

    EvtMap->clear();
}


void Gas_Scorer::DrawAll(){}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......                                                                              

void Gas_Scorer::PrintAll(){
    G4cout << " MultiFunctionalDet  " << detector->GetName() << G4endl ;
    G4cout << " PrimitiveScorer " << GetName() << G4endl               ;
    G4cout << " Number of entries " << EvtMap->entries() << G4endl     ;
}