TACTICScorer.cc 3.53 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"
5 6 7 8 9

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

Warren's avatar
Warren committed
10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
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*){

  G4double* Infos = new G4double[12];
  m_Position  = aStep->GetPreStepPoint()->GetPosition();
29

Warren's avatar
Warren committed
30 31
  Infos[0] = aStep->GetTrack()->GetTrackID();

32 33
  Infos[1] = aStep->GetTrack()->GetParticleDefinition()->GetAtomicNumber();;
  
Warren's avatar
Warren committed
34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49
  Infos[2] = aStep->GetPreStepPoint()->GetGlobalTime();
  Infos[3] = aStep->GetPreStepPoint()->GetKineticEnergy();
  Infos[4] = aStep->GetTotalEnergyDeposit() - aStep->GetNonIonizingEnergyDeposit();

  m_SegmentNumber = (int)((m_Position.z() + m_ScorerLength / 2.) / m_SegmentLength ) + 1; //Pad number 
  Infos[5] = m_SegmentNumber;
  Infos[6] = m_Position.z();
  Infos[7] = pow(pow(m_Position.x(),2) + pow(m_Position.y(),2),0.5); //R
  Infos[8] = aStep->GetTrack()->GetVertexPosition()[2];
  Infos[9] = aStep->GetTrack()->GetVertexKineticEnergy();
  G4ThreeVector p_vec = aStep->GetTrack()->GetVertexMomentumDirection();
  Infos[10] = 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[11] = aStep->GetTrack()->GetTrackLength();

  m_DetectorNumber = aStep->GetPreStepPoint()->GetTouchableHandle()->GetCopyNumber(m_Level);
  m_Index = m_DetectorNumber * 1e3 + m_SegmentNumber * 1e6;
50
  
51 52 53 54 55
  if(isnan(Infos[9])) {
    aStep->GetTrack()->SetTrackStatus(fStopAndKill);
    return 0;
  }
            
56 57 58 59 60 61
#ifdef USE_Garfield
  G4ThreeVector delta_Position = aStep->GetDeltaPosition();

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

Warren's avatar
Warren committed
72 73 74 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
  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     ;
}