diff --git a/Inputs/EventGenerator/neutron.source b/Inputs/EventGenerator/neutron.source index 4ed71397df23e884f9faac812868aab92209ed93..4f7489c76e55b60726b1f49c96fea9e973686e13 100644 --- a/Inputs/EventGenerator/neutron.source +++ b/Inputs/EventGenerator/neutron.source @@ -4,10 +4,10 @@ % Energy are given in MeV , Position in mm % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Isotropic - EnergyLow= 4 - EnergyHigh= 15 - %EnergyDistribution= flat - EnergyDistribution= Watt + EnergyLow= 0.1 + EnergyHigh= 12 + EnergyDistribution= flat + %EnergyDistribution= Watt %EnergyDistribution= 0.38*sqrt(x)*exp(-x/0.847212) %EnergyDistribution= -0.00372440431*pow(x,6)+0.387617479*pow(x,5)-14.3752948*pow(x,4)+225.888082*pow(x,3)-1555.60583*pow(x,2)+7983.24902*pow(x,1)+9069.96435 %EnergyDistribution= 0.619676*TMath::SinH(sqrt(1.07777*x))*exp(-0.847212*x) diff --git a/NPSimulation/Detectors/Vendeta/Vendeta.cc b/NPSimulation/Detectors/Vendeta/Vendeta.cc index a45b9c1f21233acc20340b43568b2ad214ab56d4..acbb987d55bce4d974a2bcb9fc1809851e8dbcbd 100644 --- a/NPSimulation/Detectors/Vendeta/Vendeta.cc +++ b/NPSimulation/Detectors/Vendeta/Vendeta.cc @@ -61,7 +61,7 @@ namespace Vendeta_NS{ // Energy and time Resolution const double EnergyThreshold = 0.1*MeV; const double ResoTime = 0.6*ns ; - const double ResoEnergy = 1.0*MeV ; + const double ResoEnergy = 0.1*MeV ; const double Thickness = 51.*mm ; const double Radius = 127./2*mm ; } @@ -345,22 +345,21 @@ void Vendeta::ReadSensitive(const G4Event* ){ double Time = RandGauss::shoot(Scorer->GetTime(i),Vendeta_NS::ResoTime); int DetectorNbr = level[0]-1; - if(Energy < 1) { - m_Event->SetHGDetectorNbr(DetectorNbr); - m_Event->SetHGQ1(0.2*Energy); - m_Event->SetHGQ2(Energy); - m_Event->SetHGTime(Time); - m_Event->SetHGQmax(0); - m_Event->SetHGIsSat(0); - } - else if(Energy > 1){ - m_Event->SetLGDetectorNbr(DetectorNbr); - m_Event->SetLGQ1(0.2*Energy); - m_Event->SetLGQ2(Energy); - m_Event->SetLGTime(Time); - m_Event->SetLGQmax(0); - m_Event->SetLGIsSat(0); - } + // Filling HG + m_Event->SetHGDetectorNbr(DetectorNbr); + m_Event->SetHGQ1(Energy); + m_Event->SetHGQ2(Energy); + m_Event->SetHGTime(Time); + m_Event->SetHGQmax(0); + m_Event->SetHGIsSat(0); + + // Filling LG + m_Event->SetLGDetectorNbr(DetectorNbr); + m_Event->SetLGQ1(Energy); + m_Event->SetLGQ2(Energy); + m_Event->SetLGTime(Time); + m_Event->SetLGQmax(0); + m_Event->SetLGIsSat(0); } } diff --git a/Projects/Vendeta_sim/Analysis.cxx b/Projects/Vendeta_sim/Analysis.cxx new file mode 100644 index 0000000000000000000000000000000000000000..dae0233c50aea05b6c406d1a48efbe85ff52b5a8 --- /dev/null +++ b/Projects/Vendeta_sim/Analysis.cxx @@ -0,0 +1,292 @@ +/***************************************************************************** + * Copyright (C) 2009-2016 this file is part of the NPTool Project * + * * + * For the licensing terms see $NPTOOL/Licence/NPTool_Licence * + * For the list of contributors see $NPTOOL/Licence/Contributors * + *****************************************************************************/ + +/***************************************************************************** + * Original Author: XAUTHORX contact address: XMAILX * + * * + * Creation Date : XMONTHX XYEARX * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * This class describe Vendeta analysis project * + * * + *---------------------------------------------------------------------------* + * Comment: * + * * + *****************************************************************************/ + +#include<iostream> +using namespace std; +#include"Analysis.h" +#include"NPAnalysisFactory.h" +#include"NPDetectorManager.h" +#include"NPOptionManager.h" +#include"RootInput.h" + +//////////////////////////////////////////////////////////////////////////////// +Analysis::Analysis(){ +} +//////////////////////////////////////////////////////////////////////////////// +Analysis::~Analysis(){ +} + +//////////////////////////////////////////////////////////////////////////////// +void Analysis::Init(){ + + Vendeta= (TVendetaPhysics*) m_DetectorManager->GetDetector("Vendeta"); + FC= (TFissionChamberPhysics*) m_DetectorManager->GetDetector("FissionChamber"); + InitialConditions = new TInitialConditions(); + + InitInputBranch(); + InitOutputBranch(); + + neutron = new NPL::Particle("1n"); +} + +//////////////////////////////////////////////////////////////////////////////// +void Analysis::TreatEvent(){ + + ReInitValue(); + + double incomingE = InitialConditions->GetKineticEnergy(0); + inEnergy = incomingE; + + unsigned int FC_mult = FC->AnodeNumber.size(); + unsigned int HF_mult = FC->Time_HF.size(); + + double GammaOffset[11]={0,0,0,0,0,0,0,0,0,0,0}; + double FakeFission_Offset = 0; + + double incomingDT=0; + double flight_path = 21500.; + if(FC_mult>0){ + int EventMax=0; + + int anode = FC->AnodeNumber[EventMax]; + double Time_FC = FC->Time[EventMax]; + bool isFake = FC->isFakeFission[EventMax]; + double Time_HF = FC->Time_HF[EventMax]; + incomingDT = FC->Time[EventMax] - FC->Time_HF[EventMax]; + TVector3 AnodePos = FC->GetVectorAnodePosition(6); + + if(anode>0){ + incomingDT -= GammaOffset[anode-1]; + AnodePos = FC->GetVectorAnodePosition(anode); + } + else if(anode ==0){ + incomingDT -= FakeFission_Offset; + } + + if(incomingDT<0){ + incomingDT += 1790; + } + + double length = flight_path;// + 6*FC->AnodeNumber[i]; + neutron->SetBeta((length/incomingDT) / c_light); + + inToF.push_back(incomingDT); + + FC_Q1.push_back(FC->Q1[EventMax]); + FC_Q2.push_back(FC->Q2[EventMax]); + FC_Qmax.push_back(FC->Qmax[EventMax]); + FC_DT.push_back(FC->DT_FC[EventMax]); + FC_FakeFission.push_back(FC->isFakeFission[EventMax]); + FC_Anode_ID.push_back(anode); + } + Vendeta->SetAnodeNumber(6); + Vendeta->BuildPhysicalEvent(); + + // VENDETA LG + unsigned int Vendeta_LG_mult = Vendeta->LG_DetectorNumber.size(); + for(unsigned int i=0; i<Vendeta_LG_mult; i++){ + + int DetNbr = Vendeta->LG_DetectorNumber[i]; + double Time_Vendeta = Vendeta->LG_Time[i]; + TVector3 DetPos = Vendeta->GetVectorDetectorPosition(DetNbr); + double Rdet = DetPos.Mag(); + //double DT = Time_Vendeta - Time_FC;// + ToF_Shift_Vendlg[DetNbr-1]; + double DT = Time_Vendeta; + + double DeltaTheta = atan(63.5/Rdet); + double Theta_Vendeta = DetPos.Theta(); + double Theta_random = ra.Uniform(Theta_Vendeta-DeltaTheta,Theta_Vendeta+DeltaTheta); + //neutron->SetTimeOfFlight(DT*1e-9/(Rdet*1e-3)); + //neutron->SetTimeOfFlight(DT*1e-9/(0.55)); + neutron->SetBeta( (Rdet/DT) / c_light); + + double En = neutron->GetEnergy(); + + // Filling output tree + LG_Tof.push_back(DT); + LG_ID.push_back(DetNbr); + LG_ELab.push_back(En); + LG_ThetaLab.push_back(Theta_random); + LG_Q1.push_back(Vendeta->LG_Q1[i]); + LG_Q2.push_back(Vendeta->LG_Q2[i]); + LG_Qmax.push_back(Vendeta->LG_Qmax[i]); + + } + + // VENDETA HG + unsigned int Vendeta_HG_mult = Vendeta->HG_DetectorNumber.size(); + for(unsigned int i=0; i<Vendeta_HG_mult; i++){ + int DetNbr = Vendeta->HG_DetectorNumber[i]; + double Time_Vendeta = Vendeta->HG_Time[i]; + TVector3 DetPos = Vendeta->GetVectorDetectorPosition(DetNbr); + double Rdet = DetPos.Mag(); + //double DT = Time_Vendeta - Time_FC;// + ToF_Shift_Vendhg[DetNbr-1]; + double DT = Time_Vendeta; + + + double DeltaTheta = atan(63.5/Rdet); + + double Theta_Vendeta = DetPos.Theta(); + double Theta_random = ra.Uniform(Theta_Vendeta-DeltaTheta,Theta_Vendeta+DeltaTheta); + //neutron->SetTimeOfFlight(DT*1e-9/(Rdet*1e-3)); + //neutron->SetTimeOfFlight(DT*1e-9/(0.55)); + neutron->SetBeta( (Rdet/DT) / c_light); + double En = neutron->GetEnergy(); + // Filling output tree + HG_ID.push_back(DetNbr); + HG_Tof.push_back(DT); + HG_ELab.push_back(En); + HG_ThetaLab.push_back(Theta_random); + HG_Q1.push_back(Vendeta->HG_Q1[i]); + HG_Q2.push_back(Vendeta->HG_Q2[i]); + HG_Qmax.push_back(Vendeta->HG_Qmax[i]); + } + + //Highlight saturated detectors + static vector<int> LG_Saturated, HG_Saturated, LG_T_sat, HG_T_sat; + LG_Saturated.clear(), HG_Saturated.clear(), LG_T_sat.clear(), HG_T_sat.clear(); + for(int j = 0; j < LG_Tof.size();j++){ + int lgID = LG_ID[j]; + if( (lgID %4 < 2) && LG_Q1[j]>500e3){ // 50 Ohm + LG_Saturated.push_back(lgID); + LG_T_sat.push_back(LG_Tof[j]); + } + else if( lgID %4 > 1 && LG_Q1[j]>590e3){ // 70 Ohm + LG_Saturated.push_back(lgID); + LG_T_sat.push_back(LG_Tof[j]); + } + } + for(int j = 0; j < HG_Tof.size();j++){ + int hgID = HG_ID[j] ; + if( (hgID %4 ==0 || hgID %4 ==1) && HG_Qmax[j] > 23800){ + HG_Saturated.push_back(hgID); + HG_T_sat.push_back(HG_Tof[j]); + } + else if( hgID %4 !=0 && hgID %4 !=1 && HG_Qmax[j] > 24300){ + HG_Saturated.push_back(hgID); + HG_T_sat.push_back(HG_Tof[j]); + } + } +} + +//////////////////////////////////////////////////////////////////////////////// +void Analysis::InitInputBranch(){ + RootInput::getInstance()->GetChain()->SetBranchStatus("InitialConditions",true); + RootInput::getInstance()->GetChain()->SetBranchStatus("fIC_*",true); + RootInput::getInstance()->GetChain()->SetBranchAddress("InitialConditions",&InitialConditions); +} + +//////////////////////////////////////////////////////////////////////////////// +void Analysis::InitOutputBranch(){ + // Incoming neutron + RootOutput::getInstance()->GetTree()->Branch("inToF",&inToF); + RootOutput::getInstance()->GetTree()->Branch("inEnergy",&inEnergy,"inEnergy/D"); + + // FissionChamber + RootOutput::getInstance()->GetTree()->Branch("FC_Q1",&FC_Q1); + RootOutput::getInstance()->GetTree()->Branch("FC_Q2",&FC_Q2); + RootOutput::getInstance()->GetTree()->Branch("FC_Qmax",&FC_Qmax); + RootOutput::getInstance()->GetTree()->Branch("FC_DT",&FC_DT); + RootOutput::getInstance()->GetTree()->Branch("FC_FakeFission",&FC_FakeFission); + RootOutput::getInstance()->GetTree()->Branch("FC_Anode_ID",&FC_Anode_ID); + /* RootOutput::getInstance()->GetTree()->Branch("FC_Time",&FC_Time); */ + /* RootOutput::getInstance()->GetTree()->Branch("HF_Time",&HF_Time); */ + + // LG + RootOutput::getInstance()->GetTree()->Branch("LG_ID",&LG_ID); + RootOutput::getInstance()->GetTree()->Branch("LG_ThetaLab",&LG_ThetaLab); + RootOutput::getInstance()->GetTree()->Branch("LG_ELab",&LG_ELab); + RootOutput::getInstance()->GetTree()->Branch("LG_Tof",&LG_Tof); + RootOutput::getInstance()->GetTree()->Branch("LG_Q1",&LG_Q1); + RootOutput::getInstance()->GetTree()->Branch("LG_Q2",&LG_Q2); + RootOutput::getInstance()->GetTree()->Branch("LG_Qmax",&LG_Qmax); + /* RootOutput::getInstance()->GetTree()->Branch("LG_Time",&LG_Time); */ + + // HG + RootOutput::getInstance()->GetTree()->Branch("HG_ID",&HG_ID); + RootOutput::getInstance()->GetTree()->Branch("HG_ThetaLab",&HG_ThetaLab); + RootOutput::getInstance()->GetTree()->Branch("HG_ELab",&HG_ELab); + RootOutput::getInstance()->GetTree()->Branch("HG_Tof",&HG_Tof); + RootOutput::getInstance()->GetTree()->Branch("HG_Q1",&HG_Q1); + RootOutput::getInstance()->GetTree()->Branch("HG_Q2",&HG_Q2); + RootOutput::getInstance()->GetTree()->Branch("HG_Qmax",&HG_Qmax); + /* RootOutput::getInstance()->GetTree()->Branch("HG_Time",&HG_Time); */ +} + +//////////////////////////////////////////////////////////////////////////////// +void Analysis::ReInitValue(){ + inToF.clear(); + inEnergy = 0; + + LG_ThetaLab.clear(); + LG_ELab.clear(); + LG_Tof.clear(); + LG_ID.clear(); + LG_Q1.clear(); + LG_Q2.clear(); + LG_Qmax.clear(); + /* LG_Time.clear(); */ + + HG_ThetaLab.clear(); + HG_ELab.clear(); + HG_Tof.clear(); + HG_ID.clear(); + HG_Q1.clear(); + HG_Q2.clear(); + HG_Qmax.clear(); + /* HG_Time.clear(); */ + + FC_Q1.clear(); + FC_Q2.clear(); + FC_Qmax.clear(); + FC_DT.clear(); + FC_FakeFission.clear(); + FC_Anode_ID.clear(); + /* FC_Time.clear(); */ + /* HF_Time.clear(); */ +} + +//////////////////////////////////////////////////////////////////////////////// +void Analysis::End(){ +} + + +//////////////////////////////////////////////////////////////////////////////// +// Construct Method to be pass to the DetectorFactory // +//////////////////////////////////////////////////////////////////////////////// +NPL::VAnalysis* Analysis::Construct(){ + return (NPL::VAnalysis*) new Analysis(); +} + +//////////////////////////////////////////////////////////////////////////////// +// Registering the construct method to the factory // +//////////////////////////////////////////////////////////////////////////////// +extern "C"{ + class proxy{ + public: + proxy(){ + NPL::AnalysisFactory::getInstance()->SetConstructor(Analysis::Construct); + } + }; + + proxy p; +} + diff --git a/Projects/Vendeta_sim/Analysis.h b/Projects/Vendeta_sim/Analysis.h new file mode 100644 index 0000000000000000000000000000000000000000..b1631630439ecbf5ad525f9334d7d6e8ef71e156 --- /dev/null +++ b/Projects/Vendeta_sim/Analysis.h @@ -0,0 +1,82 @@ +#ifndef Analysis_h +#define Analysis_h +/***************************************************************************** + * Copyright (C) 2009-2016 this file is part of the NPTool Project * + * * + * For the licensing terms see $NPTOOL/Licence/NPTool_Licence * + * For the list of contributors see $NPTOOL/Licence/Contributors * + *****************************************************************************/ + +/***************************************************************************** + * Original Author: XAUTHORX contact address: XMAILX * + * * + * Creation Date : XMONTHX XYEARX * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * This class describe Vendeta analysis project * + * * + *---------------------------------------------------------------------------* + * Comment: * + * * + *****************************************************************************/ + +#include"NPVAnalysis.h" +#include"TVendetaPhysics.h" +#include"TFissionChamberPhysics.h" +#include "NPParticle.h" +#include "TRandom3.h" +#include "TInitialConditions.h" + +class Analysis: public NPL::VAnalysis{ + public: + Analysis(); + ~Analysis(); + + public: + void Init(); + void TreatEvent(); + void End(); + void InitInputBranch(); + void InitOutputBranch(); + void ReInitValue(); + + static NPL::VAnalysis* Construct(); + + private: + vector<double> inToF; + double inEnergy; + + vector<double> LG_ID; + vector<double> LG_ThetaLab; + vector<double> LG_ELab; + vector<double> LG_Tof; + vector<double> LG_Q1; + vector<double> LG_Q2; + vector<double> LG_Qmax; + + vector<double> HG_ID; + vector<double> HG_ThetaLab; + vector<double> HG_ELab; + vector<double> HG_Tof; + vector<double> HG_Q1; + vector<double> HG_Q2; + vector<double> HG_Qmax; + + vector<double> FC_Q1; + vector<double> FC_Q2; + vector<double> FC_Qmax; + vector<double> FC_DT; + vector<bool> FC_FakeFission; + vector<int> FC_Anode_ID; + + + private: + TVendetaPhysics* Vendeta; + TFissionChamberPhysics* FC; + TInitialConditions* InitialConditions; + + NPL::Particle* neutron; + TRandom3 ra; +}; +#endif diff --git a/Projects/Vendeta_sim/CMakeLists.txt b/Projects/Vendeta_sim/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..22c74affdfc45019bdda2594f8439c52d4ab97ec --- /dev/null +++ b/Projects/Vendeta_sim/CMakeLists.txt @@ -0,0 +1,5 @@ +cmake_minimum_required (VERSION 2.8) +# Setting the policy to match Cmake version +cmake_policy(VERSION ${CMAKE_MAJOR_VERSION}.${CMAKE_MINOR_VERSION}) +# include the default NPAnalysis cmake file +include("../../NPLib/ressources/CMake/NPAnalysis.cmake") diff --git a/Projects/Vendeta_sim/PhysicsListOption.txt b/Projects/Vendeta_sim/PhysicsListOption.txt new file mode 100644 index 0000000000000000000000000000000000000000..ec1a9581f1eaf7645d75dab4f7a5ce99feb20aad --- /dev/null +++ b/Projects/Vendeta_sim/PhysicsListOption.txt @@ -0,0 +1,15 @@ +EmPhysicsList Option4 +DefaultCutOff 1 +DriftElectronPhysics 0 +IonBinaryCascadePhysics 0 +NPIonInelasticPhysics 0 +EmExtraPhysics 0 +HadronElasticPhysics 0 +StoppingPhysics 0 +OpticalPhysics 0 +HadronPhysicsINCLXX 0 +HadronPhysicsQGSP_BIC_HP 0 +HadronPhysicsQGSP_BERT_HP 0 +Decay 0 +Menate_R 0 +NeutronHP 1 diff --git a/Projects/Vendeta_sim/RunToTreat.txt b/Projects/Vendeta_sim/RunToTreat.txt new file mode 100644 index 0000000000000000000000000000000000000000..723a54e6786d8f266c9b9cd5bc0eef386e31d679 --- /dev/null +++ b/Projects/Vendeta_sim/RunToTreat.txt @@ -0,0 +1,5 @@ +TTreeName + SimulatedTree +RootFileName + ../../Outputs/Simulation/vendeta_sim_1.root + ../../Outputs/Simulation/vendeta_sim_2.root diff --git a/Projects/Vendeta_sim/Vendeta.detector b/Projects/Vendeta_sim/Vendeta.detector new file mode 100644 index 0000000000000000000000000000000000000000..8c4c5d56561d21ceb689f6011b0fc9380bccc26f --- /dev/null +++ b/Projects/Vendeta_sim/Vendeta.detector @@ -0,0 +1,18 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +Alias Theta + Action= Copy + Value= 25 36.8 48.6 60.5 72.3 84.1 95.9 107.7 119.6 131.4 143.2 155 +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +Alias Phi + Action= Copy + Value= -20 15 50 130 165 200 +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +Vendeta + R= 1000 mm + THETA= @Theta deg + PHI= @Phi deg +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +FissionChamber + POS= 0 0 0 mm + GasMaterial= CF4 + Pressure= 1 bar diff --git a/Projects/Vendeta_sim/run.mac b/Projects/Vendeta_sim/run.mac new file mode 100644 index 0000000000000000000000000000000000000000..325d88eef6b6df38e8ffac1d8fc402ccc816dc69 --- /dev/null +++ b/Projects/Vendeta_sim/run.mac @@ -0,0 +1 @@ +/run/beamOn 1000000 diff --git a/Projects/Vendeta_sim/sim.sh b/Projects/Vendeta_sim/sim.sh new file mode 100755 index 0000000000000000000000000000000000000000..054ff32b4a29d126ce0c6be8a29057193d584c7f --- /dev/null +++ b/Projects/Vendeta_sim/sim.sh @@ -0,0 +1,7 @@ +#!/bin/bash + +for i in {1..2} +do + npsimulation -D Vendeta.detector -E neutron.source --random-seed $i -O vendeta_sim_$i -B run.mac & +done +