diff --git a/Projects/Strasse/Analysis.cxx b/Projects/Strasse/Analysis.cxx index 39d571973c43584206b49115b827a2601aca6509..769562f58a39ac809d9b3a5ccb1b4163151491a5 100644 --- a/Projects/Strasse/Analysis.cxx +++ b/Projects/Strasse/Analysis.cxx @@ -45,44 +45,76 @@ void Analysis::Init(){ InitInputBranch(); Strasse = (TStrassePhysics*) m_DetectorManager -> GetDetector("Strasse"); - myReaction = new NPL::Reaction(); - myReaction->ReadConfigurationFile(NPOptionManager::getInstance()->GetReactionFile()); + // reaction properties + myQFS = new NPL::QFS(); + myQFS->ReadConfigurationFile(NPOptionManager::getInstance()->GetReactionFile()); + // reaction properties + myBeam = new NPL::Beam(); + myBeam->ReadConfigurationFile(NPOptionManager::getInstance()->GetReactionFile()); + InitialBeamEnergy = myBeam->GetEnergy() * myBeam->GetA(); // target thickness TargetThickness = m_DetectorManager->GetTargetThickness(); string TargetMaterial = m_DetectorManager->GetTargetMaterial(); + // EnergyLoss Tables + string BeamName = NPL::ChangeNameToG4Standard(myBeam->GetName()); + BeamTarget = NPL::EnergyLoss(BeamName+"_"+TargetMaterial+".G4Table","G4Table",10000); } //////////////////////////////////////////////////////////////////////////////// void Analysis::TreatEvent(){ - // Reinitiate calculated variable - ReInitValue(); - unsigned int size = Strasse->GetEventMultiplicity(); - if(size==2){ // 2 proton detected - // Proton 1 - TVector3 InnerPos1 = Strasse->GetInnerPositionOfInteraction(0); - TVector3 OuterPos1 = Strasse->GetOuterPositionOfInteraction(0); - TVector3 Proton1 = OuterPos1-InnerPos1; - // Proton 2 - TVector3 InnerPos2 = Strasse->GetInnerPositionOfInteraction(1); - TVector3 OuterPos2 = Strasse->GetOuterPositionOfInteraction(1); - TVector3 Proton2 = OuterPos2-InnerPos2; - - double deltaPhi = abs(Proton1.Phi()/deg-Proton2.Phi()/deg); - double sumTheta = Proton1.Theta()/deg+Proton2.Theta()/deg; - double OpeningAngle = Proton1.Angle(Proton2)/deg; - cout << OpeningAngle << endl; - // reject event that make no physical sense - if(deltaPhi<170 && sumTheta<80){ - return; - } - - // computing minimum distance of the two lines - TVector3 Vertex; - Distance = NPL::MinimumDistance(InnerPos1,OuterPos1,InnerPos2,OuterPos2,Vertex); - VertexX=Vertex.X(); - VertexY=Vertex.Y(); - VertexZ=Vertex.Z(); + // Reinitiate calculated variable + ReInitValue(); + unsigned int size = Strasse->GetEventMultiplicity(); + if(size==2){ // 2 proton detected + // Proton 1 + TVector3 InnerPos1 = Strasse->GetInnerPositionOfInteraction(0); + TVector3 OuterPos1 = Strasse->GetOuterPositionOfInteraction(0); + TVector3 Proton1 = OuterPos1-InnerPos1; + // Proton 2 + TVector3 InnerPos2 = Strasse->GetInnerPositionOfInteraction(1); + TVector3 OuterPos2 = Strasse->GetOuterPositionOfInteraction(1); + TVector3 Proton2 = OuterPos2-InnerPos2; + + double deltaPhi = abs(Proton1.Phi()/deg-Proton2.Phi()/deg); + double sumTheta = Proton1.Theta()/deg+Proton2.Theta()/deg; + double OpeningAngle = Proton1.Angle(Proton2)/deg; + cout << OpeningAngle << endl; + // reject event that make no physical sense + if(deltaPhi<170 && sumTheta<80){ + return; + } + + // computing minimum distance of the two lines + TVector3 Vertex; + Distance = NPL::MinimumDistance(InnerPos1,OuterPos1,InnerPos2,OuterPos2,Vertex); + VertexX=Vertex.X(); + VertexY=Vertex.Y(); + VertexZ=Vertex.Z(); } + + //double thickness_before = 0; + //double EA_vertex = BeamTarget.Slow(InitialBeamEnergy,thickness_before,0); + + // setting up Lorentz Vector from measured trajectories and energies + //LV_A.SetVect(PA); LV_p1.SetE(EA_vertex); + //LV_p1.SetVect(P1); LV_p1.SetE(E1); + //LV_p2.SetVect(P2); LV_p1.SetE(E2); + + // computing Ex from Missing Mass + //double EB = LV_A.E() + LV_T.E() - LV_p1.E() - LV_p2.E(); + //TVector3 PB = LV_A.Vect() + LV_p1.Vect() - LV_p2.Vect(); + //Ex = TMath::Sqrt( EB*EB - PB.Mag2() ) - myQFS->GetNucleusB()->Mass(); +} + + +//////////////////////////////////////////////////////////////////////////////// +TVector3 Analysis::InterpolateInPlaneZ(TVector3 V0, TVector3 V1, double Zproj){ + TVector3 Vproj(-999,-999,-999); + double t = (Zproj - V1.Z()) / (V1.Z()-V0.Z()); + double Xproj= V1.X() + (V1.X()-V0.X()) * t; + double Yproj= V1.Y() + (V1.Y()-V0.Y()) * t; + Vproj.SetXYZ(Xproj,Yproj,Zproj); + return Vproj; } //////////////////////////////////////////////////////////////////////////////// diff --git a/Projects/Strasse/strasse_optimized.detector b/Projects/Strasse/strasse_optimized.detector index 201a72f39155b875102f8d149c7404591761ed5e..4b341df72cedd879576b644698741012044c85a2 100644 --- a/Projects/Strasse/strasse_optimized.detector +++ b/Projects/Strasse/strasse_optimized.detector @@ -90,7 +90,7 @@ Strasse Inner %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%1 Strasse Outer Radius= @OuterRadius mm - Z= 91.0 mm + Z= 78.0 mm Phi= @OuterPhi deg Shift= @OuterShift mm diff --git a/Projects/Strasse/strasse_ref.detector b/Projects/Strasse/strasse_ref.detector index 8789c5858e4bc6c4e565fe2ee3046981ec74bd97..9ea0d4fc007d1ab5886a9eda3fb4e32a3cd79112 100644 --- a/Projects/Strasse/strasse_ref.detector +++ b/Projects/Strasse/strasse_ref.detector @@ -43,9 +43,17 @@ Strasse Info Outer_PCB_Thickness= 1.6 mm Outer_Wafer_TransverseStrips= 605 Outer_Wafer_LongitudinalStrips= 325 + % all radius are external, internal deduced from thickness Chamber_Thickness= 3 mm - Chamber_Length= 500 mm - Chamber_Radius= 100 mm + Chamber_Cylinder_Length= 360 mm + Chamber_Radius= 180 mm + Chamber_ExitTube_Radius= 79.5 mm + Chamber_ExitTube_Length= 100 mm + Chamber_Flange_Inner_Radius= 50 mm + Chamber_Sphere_Radius= 220 mm + Chamber_Sphere_Shift= 60 mm + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%1 Alias InnerShift @@ -83,7 +91,7 @@ Strasse Inner %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%1 Strasse Outer Radius= @OuterRadius mm - Z= 91.0 mm + Z= 78.0 mm Phi= @OuterPhi deg Shift= @OuterShift mm @@ -92,6 +100,6 @@ Strasse Chamber Z= -30 mm %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -Catana Dummy - Z= 300 mm +%Catana Dummy +% Z= 300 mm