* Prepare analysis for missing mass + update geometries

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;
 }
 
 ////////////////////////////////////////////////////////////////////////////////

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
 

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