From 3a883ca7debe19546fd06f7b3fa8935324d17d09 Mon Sep 17 00:00:00 2001 From: Warren Lynch <warren.lynch@york.ac.uk> Date: Wed, 24 Mar 2021 14:09:09 +0000 Subject: [PATCH] Please enter the commit message for your changes. Lines starting --- NPSimulation/Detectors/TACTIC/GARFDRIFT.h | 43 +++++++++++++------ NPSimulation/Detectors/TACTIC/TACTIC.cc | 7 +-- NPSimulation/Detectors/TACTIC/TACTICScorer.cc | 28 ++++++++---- 3 files changed, 53 insertions(+), 25 deletions(-) diff --git a/NPSimulation/Detectors/TACTIC/GARFDRIFT.h b/NPSimulation/Detectors/TACTIC/GARFDRIFT.h index f86771dc1..a3141e8aa 100644 --- a/NPSimulation/Detectors/TACTIC/GARFDRIFT.h +++ b/NPSimulation/Detectors/TACTIC/GARFDRIFT.h @@ -8,28 +8,34 @@ #include "Garfield/ViewSignal.hh" #include "Garfield/ViewMedium.hh" -void GARFDRIFT(double energy, double t_start, G4ThreeVector start, G4ThreeVector delta_pos, double R, int Pad_start, double ID, - double ScoLength, double SegLength) { +double GARFDRIFT(double energy, double t_start, G4ThreeVector start, G4ThreeVector delta_pos, double R, int Pad_start, double ID, + double ScoLength, double SegLength, double event) { const double rWire = 1.2; const double rTube = 5.; const double lTube = 25.19; + //TH1F *hist = new TH1F("","",1000,0,10000); // 10 ns bins //For R6 = 2, R7 = 1 //const double vWire = -375.; //For V_total = 750V //const double vWire = -645.3312; //For V_total = 1200V //SHOULD BE 583.741 !? //const double vWire = -535.0959; //For V_total = 1100V //const double vWire = -875.6116; //For V_total = 1800V - const double vWire = -681.0312; //For V_total = 1400V + //const double vWire = -681.0312; //For V_total = 1400V + //const double vWire = -730.; + const double vWire = -1000.; const double vTube = 0.; double x_start, y_start, z_start; ofstream file; - double production_bias = 1.e-01; + double production_bias = 0.001; double t_end, x_end, y_end, z_end; double Pad_end; //double w_value = 35.; //for HeCO2 --guess - double w_value = 26.31; //for argon + //double w_value = 26.31; //for argon + double w_value = 41.1;// for He int electrons = (int)(energy/w_value*production_bias);//*production_bias; //Need to carry excess energy over later. + double energy_excess = 0.; + vector<double> pad_vec, time_vec; Garfield::MediumMagboltz* gas = new Garfield::MediumMagboltz(); Garfield::ViewMedium* mediumView = new Garfield::ViewMedium(); @@ -42,10 +48,15 @@ void GARFDRIFT(double energy, double t_start, G4ThreeVector start, G4ThreeVector cout << "\n" << electrons << "\n" << endl; if(R>rWire && electrons>0) { - - gas->LoadGasFile("ar_90_ch4_10.gas"); //atomic fraction in this case (P10). + + energy_excess = (energy/w_value*production_bias - electrons)*w_value/production_bias; + + //gas->LoadGasFile("he_90_co2_10_200mbar.gas"); + //gas->LoadGasFile("ar_90_ch4_10.gas"); //atomic fraction in this case (P10). //gas->LoadGasFile("he_90_co2_10_500mbar.gas"); //mass fraction in this case - + //gas->LoadGasFile("he_90_co2_10_1bar.gas"); + gas->LoadGasFile("he_90_co2_10_100mbar.gas"); + gas->Initialise(true); mediumView->SetMedium(gas); @@ -69,14 +80,18 @@ void GARFDRIFT(double energy, double t_start, G4ThreeVector start, G4ThreeVector drift->GetDriftLinePoint(np-1, x_end, y_end, z_end, t_end); //np-1 is last DriftLineEndPoint Pad_end = (int)((z_end + ScoLength / 2.) / SegLength ) + 1; //new Pad number - - file.open("signal.dat", std::ios::app); - file << Pad_end << "\t" << t_end << endl; - file.close(); - + + if(pow(pow(x_end,2)+pow(y_end,2),0.5)>4.9) { //reached the anode + file.open("signal.dat", std::ios::app); + file << event << "\t" << ID << "\t" << Pad_start << "\t" << Pad_end << "\t" << t_end << endl; + file.close(); + } + } } - + delete gas; delete mediumView; delete geo; delete tube; delete cmp; delete sensor; delete drift; //otherwise these are held in memory and cause a killed: 9 crash. + + return energy_excess; } diff --git a/NPSimulation/Detectors/TACTIC/TACTIC.cc b/NPSimulation/Detectors/TACTIC/TACTIC.cc index 5a2989c81..c39020979 100644 --- a/NPSimulation/Detectors/TACTIC/TACTIC.cc +++ b/NPSimulation/Detectors/TACTIC/TACTIC.cc @@ -404,12 +404,13 @@ void TACTIC::InitializeScorers() { */ LightFilter->addIon(m_Reaction.GetParticle3()->GetZ(),m_Reaction.GetParticle3()->GetA()); HeavyFilter->addIon(m_Reaction.GetParticle4()->GetZ(),m_Reaction.GetParticle4()->GetA()); - LightScorer->SetFilter(LightFilter); - HeavyScorer->SetFilter(HeavyFilter); if(m_Reaction.GetParticle1()->GetZ() == m_Reaction.GetParticle4()->GetZ()) BeamFilter->add("geantino"); else BeamFilter->addIon(m_Reaction.GetParticle1()->GetZ(),m_Reaction.GetParticle1()->GetA()); - BeamScorer->SetFilter(BeamFilter); + LightScorer->SetFilter(LightFilter); + HeavyScorer->SetFilter(HeavyFilter); + BeamScorer->SetFilter(BeamFilter); + m_Scorer->RegisterPrimitive(LightScorer); m_Scorer->RegisterPrimitive(HeavyScorer); m_Scorer->RegisterPrimitive(BeamScorer); diff --git a/NPSimulation/Detectors/TACTIC/TACTICScorer.cc b/NPSimulation/Detectors/TACTIC/TACTICScorer.cc index 79e20f0d6..1bc7904db 100644 --- a/NPSimulation/Detectors/TACTIC/TACTICScorer.cc +++ b/NPSimulation/Detectors/TACTIC/TACTICScorer.cc @@ -29,7 +29,7 @@ Gas_Scorer::~Gas_Scorer(){} G4bool Gas_Scorer::ProcessHits(G4Step* aStep, G4TouchableHistory*){ - G4double* Infos = new G4double[14]; + G4double* Infos = new G4double[15]; m_Position = aStep->GetPreStepPoint()->GetPosition(); Infos[0] = G4RunManager::GetRunManager()->GetCurrentEvent()->GetEventID(); @@ -60,18 +60,30 @@ G4bool Gas_Scorer::ProcessHits(G4Step* aStep, G4TouchableHistory*){ aStep->GetTrack()->SetTrackStatus(fStopAndKill); return 0; } - -#ifdef USE_Garfield - G4ThreeVector delta_Position = aStep->GetDeltaPosition(); - - 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]); -#endif - + map<G4int, G4double**>::iterator it; it= EvtMap->GetMap()->find(m_Index); if(it!=EvtMap->GetMap()->end()){ G4double* dummy = *(it->second); if(Infos[1]==dummy[1]) Infos[5]+=dummy[5]; //accumulate ionisation energy deposit to get total accross pad + + +#ifdef USE_Garfield + G4ThreeVector delta_Position = aStep->GetDeltaPosition(); + G4double excess; + + if(aStep->IsFirstStepInVolume() == 1) excess = 0.; + if(aStep->IsFirstStepInVolume() == 0) excess = dummy[14]/eV; + + if(excess > 1.e06) excess = 0.; //If dummy[12] is not defined returns random value > 1.e06 ? + + if(Infos[8] > 12.) Infos[14] = GARFDRIFT((Infos[5]/eV+excess), Infos[3], m_Position/cm, delta_Position/cm, Infos[8]/cm, Infos[6], Infos[2], m_ScorerLength/cm, m_SegmentLength/cm, Infos[0])*eV; +#endif + if(Infos[8] < 12.) Infos[14] = 0; + + + if(Infos[14]>0.) cout << "Infos[14] " << Infos[14]/eV << endl; + delete dummy; } -- GitLab