From 2d388af87e98176f92dda0e764061b34ea358433 Mon Sep 17 00:00:00 2001
From: adrien matta <matta@lpccaen.in2p3.fr>
Date: Mon, 20 Jul 2020 16:48:38 +0200
Subject: [PATCH] * Strasse position now working !

---
 NPLib/Detectors/Strasse/TStrassePhysics.cxx | 112 ++++++++++++--------
 NPSimulation/Detectors/Strasse/Strasse.cc   |   8 +-
 Projects/Strasse/Analysis.cxx               |   4 +
 Projects/Strasse/Analysis.h                 |   2 +
 Projects/Strasse/Control.cxx                |  33 ++++++
 Projects/Strasse/PhysicsListOption.txt      |   2 +-
 6 files changed, 110 insertions(+), 51 deletions(-)
 create mode 100644 Projects/Strasse/Control.cxx

diff --git a/NPLib/Detectors/Strasse/TStrassePhysics.cxx b/NPLib/Detectors/Strasse/TStrassePhysics.cxx
index ea2aff7d3..54fb6719c 100644
--- a/NPLib/Detectors/Strasse/TStrassePhysics.cxx
+++ b/NPLib/Detectors/Strasse/TStrassePhysics.cxx
@@ -55,7 +55,7 @@ ClassImp(TStrassePhysics)
     m_NumberOfInnerDetectors = 0;
     m_NumberOfOuterDetectors = 0;
     m_MaximumStripMultiplicityAllowed = 10;
-    m_StripEnergyMatching = 0.050;
+    m_StripEnergyMatching = 1.00;
 
     ////////////////////
     // Inner Detector //
@@ -110,18 +110,18 @@ ClassImp(TStrassePhysics)
 void TStrassePhysics::AddInnerDetector(double R, double Z, double Phi, double Shift){
   m_NumberOfInnerDetectors++;
   double ActiveWidth  = Inner_Wafer_Width-2.*Inner_Wafer_GuardRing;
-  double ActiveLength = Inner_Wafer_Length-Inner_Wafer_PADExternal-Inner_Wafer_PADInternal-Inner_Wafer_GuardRing;
-  double LongitudinalPitch = ActiveLength/Inner_Wafer_TransverseStrips;
-  double TransversePitch = ActiveWidth/Inner_Wafer_LongitudinalStrips;
-  TVector3 Det_pos(Shift,R,Z) ;
+  double ActiveLength = Inner_Wafer_Length-Inner_Wafer_PADExternal-Inner_Wafer_PADInternal-2*Inner_Wafer_GuardRing;
+  double LongitudinalPitch = ActiveWidth/Inner_Wafer_LongitudinalStrips;
+  double TransversePitch = ActiveLength/Inner_Wafer_TransverseStrips;
+//cout << ActiveWidth << " " << ActiveLength << " " << LongitudinalPitch << " " << TransversePitch << endl;
 
   // Vector C position of detector face center
   TVector3 C(Shift,R,Z);// center of the whole detector, including PCB
-  C.RotateZ(Phi);
+  C.RotateZ(-Phi);
 
   // Vector W normal to detector face (pointing to the back)
   TVector3 W(0,1,0);
-  W.RotateZ(Phi);
+  W.RotateZ(-Phi);
 
   // Vector U on detector face (parallel to Z axis/longitudinal strips) 
   TVector3 U = TVector3(0,0,1);
@@ -131,11 +131,16 @@ void TStrassePhysics::AddInnerDetector(double R, double Z, double Phi, double Sh
   // Adding position for downstream silicon:
   // Moving to corner of the silicon
   TVector3 P_1_1 = C
-    +U*(Inner_PCB_UpstreamWidth-Inner_PCB_DownstreamWidth)
-    +U*(0.5*Inner_PCB_MidWidth+Inner_Wafer_GuardRing+Inner_Wafer_PADInternal)
-    +V*(Inner_PCB_PortWidth-Inner_PCB_StarboardWidth)
-    +V*(0.5*Inner_Wafer_Width-Inner_Wafer_GuardRing);
-
+      +U*0.5*(Inner_PCB_UpstreamWidth-Inner_PCB_DownstreamWidth) // In between wafer
+      -U*0.5*Inner_PCB_MidWidth // Internal wafer edge
+      -U*Inner_Wafer_Length // External wafer edge
+      +U*(Inner_Wafer_GuardRing+Inner_Wafer_PADExternal) // External active wafer edge
+      +U*0.5*TransversePitch // middle of strip
+      -V*0.5*(Inner_PCB_StarboardWidth-Inner_PCB_PortWidth)
+      -V*0.5*Inner_Wafer_Width
+      +V*Inner_Wafer_GuardRing
+      +V*0.5*LongitudinalPitch; // middle of strip
+  
   vector<double> lineX;
   vector<double> lineY;
   vector<double> lineZ;
@@ -146,19 +151,19 @@ void TStrassePhysics::AddInnerDetector(double R, double Z, double Phi, double Sh
 
   TVector3 P;
   for(int i=0; i<Inner_Wafer_TransverseStrips; i++){
-  lineX.clear();
-  lineY.clear();
-  lineZ.clear();
-  for(int j=0; j<Inner_Wafer_LongitudinalStrips; j++){
-  P = P_1_1 + i*U*LongitudinalPitch + j*V*TransversePitch;
-  lineX.push_back(P.X());
-  lineY.push_back(P.Y());
-  lineZ.push_back(P.Z());
-  }
+    lineX.clear();
+    lineY.clear();
+    lineZ.clear();
+    for(int j=0; j<Inner_Wafer_LongitudinalStrips; j++){
+      P = P_1_1 + i*U*TransversePitch + j*V*LongitudinalPitch;
+      lineX.push_back(P.X());
+      lineY.push_back(P.Y());
+      lineZ.push_back(P.Z());
+      }
 
-  OneDetectorStripPositionX.push_back(lineX);
-  OneDetectorStripPositionY.push_back(lineY);
-  OneDetectorStripPositionZ.push_back(lineZ);
+    OneDetectorStripPositionX.push_back(lineX);
+    OneDetectorStripPositionY.push_back(lineY);
+    OneDetectorStripPositionZ.push_back(lineZ);
   }
 
   m_InnerStripPositionX.push_back(OneDetectorStripPositionX);
@@ -168,18 +173,32 @@ void TStrassePhysics::AddInnerDetector(double R, double Z, double Phi, double Sh
   // Adding position for upstream silicon:
   // Moving to corner of the silicon
   P_1_1 = C
-    -U*(Inner_PCB_UpstreamWidth-Inner_PCB_DownstreamWidth)
-    -U*(0.5*Inner_PCB_MidWidth+Inner_Wafer_GuardRing+Inner_Wafer_PADInternal)
-    +V*(Inner_PCB_PortWidth-Inner_PCB_StarboardWidth)
-    +V*(0.5*Inner_Wafer_Width-Inner_Wafer_GuardRing);
+      +U*0.5*(Inner_PCB_UpstreamWidth-Inner_PCB_DownstreamWidth) // In between wafer
+      +U*0.5*Inner_PCB_MidWidth // Internal wafer edge
+      +U*(Inner_Wafer_GuardRing+Inner_Wafer_PADInternal) // Internal active wafer edge
+      +U*0.5*TransversePitch// middle of strip
+      -V*0.5*(Inner_PCB_StarboardWidth-Inner_PCB_PortWidth)
+      -V*0.5*Inner_Wafer_Width
+      +V*Inner_Wafer_GuardRing
+      +V*0.5*LongitudinalPitch; // middle of strip
 
   for(int i=0; i<Inner_Wafer_TransverseStrips; i++){
-  for(int j=0; j<Inner_Wafer_LongitudinalStrips; j++){
-    P = P_1_1 + i*U*LongitudinalPitch - j*V*TransversePitch;
-    m_InnerStripPositionX[m_NumberOfInnerDetectors-1][i].push_back(P.X());
-    m_InnerStripPositionY[m_NumberOfInnerDetectors-1][i].push_back(P.Y());
-    m_InnerStripPositionZ[m_NumberOfInnerDetectors-1][i].push_back(P.Z());
+    lineX.clear();
+    lineY.clear();
+    lineZ.clear();
+
+    for(int j=0; j<Inner_Wafer_LongitudinalStrips; j++){
+      P = P_1_1 + i*U*TransversePitch + j*V*LongitudinalPitch;
+      lineX.push_back(P.X());
+      lineY.push_back(P.Y());
+      lineZ.push_back(P.Z());
+
     }
+
+    m_InnerStripPositionX[m_NumberOfInnerDetectors-1].push_back(lineX);
+    m_InnerStripPositionY[m_NumberOfInnerDetectors-1].push_back(lineY);
+    m_InnerStripPositionZ[m_NumberOfInnerDetectors-1].push_back(lineZ);
+
   }
 }
 
@@ -189,9 +208,10 @@ void TStrassePhysics::AddOuterDetector(double R, double Z, double Phi, double Sh
 
 ///////////////////////////////////////////////////////////////////////////
 TVector3 TStrassePhysics::GetInnerPositionOfInteraction(const int i){
-  TVector3 Position = TVector3(GetInnerStripPositionX(DetectorNumber[i], InnerStripL[i], InnerStripT[i]),
-      GetInnerStripPositionY(DetectorNumber[i], InnerStripL[i], InnerStripT[i]),
-      GetInnerStripPositionZ(DetectorNumber[i], InnerStripL[i], InnerStripT[i]));
+  TVector3 Position = TVector3(
+      GetInnerStripPositionX(DetectorNumber[i], InnerStripT[i], InnerStripL[i]),
+      GetInnerStripPositionY(DetectorNumber[i], InnerStripT[i], InnerStripL[i]),
+      GetInnerStripPositionZ(DetectorNumber[i], InnerStripT[i], InnerStripL[i]));
 
   return Position;
 }
@@ -289,8 +309,6 @@ vector<TVector2> TStrassePhysics::Match_X_Y(){
         // Declaration of variable for clarity
         double TE = m_PreTreatedData->GetInner_TE_Energy(i);
         double LE = m_PreTreatedData->GetInner_LE_Energy(i);
-        double XStripNbr = m_PreTreatedData->GetInner_TE_StripNbr(i);
-        double YStripNbr = m_PreTreatedData->GetInner_LE_StripNbr(i);
 
         // look if energy matches
         if(abs(TE-LE)/2.<m_StripEnergyMatching){
@@ -443,14 +461,6 @@ void TStrassePhysics::ReadAnalysisConfig() {
 void TStrassePhysics::Clear() {
   EventMultiplicity = 0;
 
-  // Position Information
-  InnerPosX.clear();
-  InnerPosY.clear();
-  InnerPosZ.clear();
-  OuterPosX.clear();
-  OuterPosY.clear();
-  OuterPosZ.clear();
-
   // DSSD
   DetectorNumber.clear();
   E.clear();
@@ -459,6 +469,16 @@ void TStrassePhysics::Clear() {
   OuterStripT.clear();
   OuterStripL.clear();
   DE.clear();
+
+  // Position Information
+  InnerPosX.clear();
+  InnerPosY.clear();
+  InnerPosZ.clear();
+  OuterPosX.clear();
+  OuterPosY.clear();
+  OuterPosZ.clear();
+
+
 }
 
 
diff --git a/NPSimulation/Detectors/Strasse/Strasse.cc b/NPSimulation/Detectors/Strasse/Strasse.cc
index 292ff57a5..24f745fc2 100644
--- a/NPSimulation/Detectors/Strasse/Strasse.cc
+++ b/NPSimulation/Detectors/Strasse/Strasse.cc
@@ -585,7 +585,7 @@ void Strasse::ConstructDetector(G4LogicalVolume* world){
  
   // Inner Barrel
   for (unsigned short i = 0 ; i < m_Inner_R.size() ; i++) {
-    G4ThreeVector Det_pos = G4ThreeVector(m_Inner_Shift[i],m_Inner_R[i],m_Inner_Z[i]) ;
+    G4ThreeVector Det_pos = G4ThreeVector(m_Inner_Shift[i],m_Inner_R[i]+0.5*Inner_PCB_Thickness,m_Inner_Z[i]) ;
     Det_pos.rotate(-m_Inner_Phi[i],G4ThreeVector(0,0,1));
     G4RotationMatrix* Rot =  new G4RotationMatrix(0*deg,0*deg,m_Inner_Phi[i]);
 
@@ -596,7 +596,7 @@ void Strasse::ConstructDetector(G4LogicalVolume* world){
 
   // Outer Barrel 
   for (unsigned short i = 0 ; i < m_Outer_R.size() ; i++) {
-    G4ThreeVector Det_pos = G4ThreeVector(m_Outer_Shift[i],m_Outer_R[i],m_Outer_Z[i]) ;
+    G4ThreeVector Det_pos = G4ThreeVector(m_Outer_Shift[i],m_Outer_R[i]+0.5*Inner_PCB_Thickness,m_Outer_Z[i]) ;
     Det_pos.rotate(-m_Outer_Phi[i],G4ThreeVector(0,0,1));
     G4RotationMatrix* Rot =  new G4RotationMatrix(0*deg,0*deg,m_Outer_Phi[i]);
 
@@ -741,7 +741,7 @@ void Strasse::InitializeScorers() {
   // Otherwise the scorer is initialised
   m_Active_InnerWafer_Width= Inner_Wafer_Width-2.*Inner_Wafer_GuardRing;
   m_Active_InnerWafer_Length= 
-  Inner_Wafer_Length-Inner_Wafer_PADExternal-Inner_Wafer_PADInternal-Inner_Wafer_GuardRing;
+  Inner_Wafer_Length-Inner_Wafer_PADExternal-Inner_Wafer_PADInternal-2*Inner_Wafer_GuardRing;
     
 
   G4VPrimitiveScorer* InnerScorer1 = new DSSDScorers::PS_Rectangle("InnerScorer1",2,
@@ -760,7 +760,7 @@ void Strasse::InitializeScorers() {
 
   m_Active_OuterWafer_Width=Outer_Wafer_Width-2.*Outer_Wafer_GuardRing;
   m_Active_OuterWafer_Length=
-  Outer_Wafer_Length-Outer_Wafer_PADExternal-Outer_Wafer_PADInternal-Outer_Wafer_GuardRing;
+  Outer_Wafer_Length-Outer_Wafer_PADExternal-Outer_Wafer_PADInternal-2*Outer_Wafer_GuardRing;
 
 
   G4VPrimitiveScorer* OuterScorer1 = new DSSDScorers::PS_Rectangle("OuterScorer1",2,
diff --git a/Projects/Strasse/Analysis.cxx b/Projects/Strasse/Analysis.cxx
index 115d337aa..ac3ada148 100644
--- a/Projects/Strasse/Analysis.cxx
+++ b/Projects/Strasse/Analysis.cxx
@@ -35,6 +35,7 @@ Analysis::~Analysis(){
 
 ////////////////////////////////////////////////////////////////////////////////
 void Analysis::Init(){
+  IC= new TInteractionCoordinates;
   InitOutputBranch();
   InitInputBranch();
   
@@ -61,10 +62,13 @@ void Analysis::InitOutputBranch() {
   RootOutput::getInstance()->GetTree()->Branch("ELab",&ELab,"ELab/D");
   RootOutput::getInstance()->GetTree()->Branch("ThetaLab",&ThetaLab,"ThetaLab/D");
   RootOutput::getInstance()->GetTree()->Branch("ThetaCM",&ThetaCM,"ThetaCM/D");
+  RootOutput::getInstance()->GetTree()->Branch("InteractionCoordinates","TInteractionCoordinates",&IC);
 }
 
 ////////////////////////////////////////////////////////////////////////////////
 void Analysis::InitInputBranch(){
+
+    RootInput:: getInstance()->GetChain()->SetBranchAddress("InteractionCoordinates",&IC);
 }
 ////////////////////////////////////////////////////////////////////////////////
 void Analysis::ReInitValue(){
diff --git a/Projects/Strasse/Analysis.h b/Projects/Strasse/Analysis.h
index 9492dd1f2..ff841d675 100644
--- a/Projects/Strasse/Analysis.h
+++ b/Projects/Strasse/Analysis.h
@@ -28,6 +28,7 @@
 #include"RootInput.h"
 #include "TStrassePhysics.h"
 #include "TInitialConditions.h"
+#include "TInteractionCoordinates.h"
 #include <TRandom3.h>
 #include <TVector3.h>
 #include <TMath.h>
@@ -83,5 +84,6 @@ TInitialConditions* myInit ;
   double Si_X_Sharc ;
   double Si_Y_Sharc ;
   TStrassePhysics* Strasse;
+  TInteractionCoordinates* IC;
 };
 #endif
diff --git a/Projects/Strasse/Control.cxx b/Projects/Strasse/Control.cxx
new file mode 100644
index 000000000..2a5a4e1eb
--- /dev/null
+++ b/Projects/Strasse/Control.cxx
@@ -0,0 +1,33 @@
+void Control(){
+  TFile* file = TFile::Open("../../Outputs/Analysis/PhysicsTree.root");
+  TTree* PhysicsTree= (TTree*) file->FindObjectAny("PhysicsTree");
+  TCanvas* c = new TCanvas("Control","Control",1000,1000);
+  c->Divide(2,2);
+  //string cond = "InnerPosX!=-60 && InnerStripL[0]==64  && InnerStripT[0]==152";
+
+  string cond = "InnerPosX!=-60 && (DetectorNumber==3 || DetectorNumber==1)";
+  c->cd(1);
+//  PhysicsTree->Draw("InnerPosX:fDetected_Position_X",cond.c_str(),"col") ; 
+//
+  PhysicsTree->Draw("InnerPosX-fDetected_Position_X>>h(1000,-2,2)",cond.c_str(),"col") ; 
+  c->cd(2);
+  //PhysicsTree->Draw("InnerPosY:fDetected_Position_Y",cond.c_str(),"col") ; 
+
+  PhysicsTree->Draw("InnerPosY-fDetected_Position_Y>>h2(1000,-2,2)",cond.c_str(),"col") ; 
+  c->cd(3);
+  //PhysicsTree->Draw("InnerPosZ:fDetected_Position_Z",cond.c_str(),"col") ; 
+
+  PhysicsTree->Draw("InnerPosZ-fDetected_Position_Z>>h3(1000,-2,2)",cond.c_str(),"col") ; 
+  c->cd(3);
+  c->cd(4);
+  PhysicsTree->Draw("InnerPosY:InnerPosX:InnerPosZ", cond.c_str(),"");            
+
+ /* TCanvas* c2 = new TCanvas("Control 2", "Control2",500,500,2000,1000);  
+  c2->Divide(2,1);
+  c2->cd(1);
+  PhysicsTree->Draw("InnerPosY:InnerPosX","InnerPosX!=-60");            
+
+  c2->cd(2);
+  PhysicsTree->Draw("fDetected_Position_Y:fDetected_Position_X","InnerPosX!=-60","");            
+*/
+}
diff --git a/Projects/Strasse/PhysicsListOption.txt b/Projects/Strasse/PhysicsListOption.txt
index 444fc7cce..ff38d11c7 100644
--- a/Projects/Strasse/PhysicsListOption.txt
+++ b/Projects/Strasse/PhysicsListOption.txt
@@ -1,5 +1,5 @@
 EmPhysicsList Option4
-DefaultCutOff 1
+DefaultCutOff 1000000
 IonBinaryCascadePhysics 0
 NPIonInelasticPhysics 0
 EmExtraPhysics 0
-- 
GitLab