diff --git a/Inputs/DetectorConfiguration/Riken_65mm.detector b/Inputs/DetectorConfiguration/Riken_65mm.detector
index 3754de046edaeb0802e0f2837bf2a5c88c56d340..efc9692767cd8e6c6c691e153d9615d0757508b1 100644
--- a/Inputs/DetectorConfiguration/Riken_65mm.detector
+++ b/Inputs/DetectorConfiguration/Riken_65mm.detector
@@ -7,7 +7,7 @@
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
GeneralTarget
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-Target
+%Target
THICKNESS= 1
RADIUS= 45
MATERIAL= CD2
@@ -83,6 +83,27 @@ SI= 1
SILI= 0
CSI= 1
VIS= all
+%%%%%%% Telescope 6 %%%%%%%
+M2Telescope
+X1_Y1= 155.77 50.23 8.18
+X1_Y128= 155.77 -50.23 8.18
+X128_Y1= 133.17 50.23 -89.7
+X128_Y128= 133.17 -50.23 -89.7
+SI= 1
+SILI= 1
+CSI= 0
+VIS= all
+
+%%%%%%% Telescope 7 %%%%%%%
+M2Telescope
+X1_Y1= 27.07 50.23 -154.49
+X1_Y128= 116.58 50.23 -108.88
+X128_Y1= 27.07 -50.23 -154.49
+X128_Y128= 116.58 -50.23 -108.88
+SI= 1
+SILI= 1
+CSI= 0
+VIS= all
%%%%%%%%%%%%%%%%%%%%
AddThinSi
diff --git a/NPAnalysis/10He_Riken/Analysis b/NPAnalysis/10He_Riken/Analysis
index a183cdf33d9798cb81ea8c62a0a7ba06212add83..6c0904b1217ac65c52f5f7efcf964871d826d228 100755
Binary files a/NPAnalysis/10He_Riken/Analysis and b/NPAnalysis/10He_Riken/Analysis differ
diff --git a/NPAnalysis/10He_Riken/include/ObjectManager.hh b/NPAnalysis/10He_Riken/include/ObjectManager.hh
index 3e2c4435f74626023ca02557ed06719e8ef50154..8d2739151bf089ed4af970bb3bb3d7f29cc583e0 100644
--- a/NPAnalysis/10He_Riken/include/ObjectManager.hh
+++ b/NPAnalysis/10He_Riken/include/ObjectManager.hh
@@ -104,7 +104,7 @@ using namespace NPL ;
namespace ENERGYLOSS
{
- // Declare your Energy loss here :
+ // 3He Energy Loss
EnergyLoss He3TargetWind = EnergyLoss ( "3He_Mylar.txt" ,
1000 ,
1 ,
@@ -124,6 +124,24 @@ namespace ENERGYLOSS
10 ,
1 ,
3 );
+
+
+ // proton Energy Loss
+ EnergyLoss protonTargetWind = EnergyLoss ( "proton_Mylar.txt" ,
+ 1000 ,
+ 1 ,
+ 1 );
+
+ EnergyLoss protonTargetGaz = EnergyLoss ( "proton_D2gaz_1b_26K.txt" ,
+ 1000 ,
+ 1 ,
+ 1 );
+
+ EnergyLoss protonStripAl = EnergyLoss ( "proton_Al.txt" ,
+ 10 ,
+ 1 ,
+ 1 );
+
}
diff --git a/NPAnalysis/10He_Riken/src/Analysis.cc b/NPAnalysis/10He_Riken/src/Analysis.cc
index 6df2fcbf98e8a79b1a4657c2d562f84483846fb4..a09db085a49b944b53913ca04d20a004871cd329 100644
--- a/NPAnalysis/10He_Riken/src/Analysis.cc
+++ b/NPAnalysis/10He_Riken/src/Analysis.cc
@@ -42,7 +42,7 @@ int main(int argc,char** argv)
RootOutput::getInstance()->GetTree()->Branch("ThetaCM",&ThetaCM,"ThetaCM/D") ;
RootOutput::getInstance()->GetTree()->Branch("ResolThetaCM",&ResolThetaCM,"ResolThetaCM/D") ;
// Get the formed Chained Tree and Treat it
- //TChain* Chain = RootInput:: getInstance() -> GetChain() ;
+ TChain* Chain = RootInput:: getInstance() -> GetChain() ;
// Open the ThinSi Branch
Chain->SetBranchStatus("ThinSiEnergy",true) ;
@@ -86,71 +86,109 @@ int main(int argc,char** argv)
double ThetaN = ThetaCalculation ( HitDirection , TVector3(0,0,1) ) ;
// ANgle between particule and Must2 Si surface
double ThetaMM2Surface = ThetaCalculation ( HitDirection , M2 -> GetTelescopeNormal() );
-
- if(E>-1000 && ThinSi>0 )
+ if(M2 -> GetPositionOfInteraction().Z()>0)
{
- E= He3StripAl.EvaluateInitialEnergy( E , // Energy of the detected particle
- 2*0.4*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
-
- E= He3StripSi.EvaluateInitialEnergy( E , // Energy of the detected particle
- 20*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
-
-// E = E + ThinSi ;
+ if(E>-1000 && ThinSi>0 )
+ {
+ E= He3StripAl.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 2*0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
+
+ E= He3StripSi.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 20*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
- E= He3StripAl.EvaluateInitialEnergy( E , // Energy of the detected particle
- 0.4*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
+ // E = E + ThinSi ;
+
+ E= He3StripAl.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
- // cout << E << " " << Eb-E << " " << ThinSi << endl ;
+ // cout << E << " " << Eb-E << " " << ThinSi << endl ;
- E= He3TargetWind.EvaluateInitialEnergy( E , // Energy of the detected particle
- 15*micrometer , // Target Thickness at 0 degree
- ThetaN );
-
- E= He3TargetGaz.EvaluateInitialEnergy( E , // Energy of the detected particle
- 1.5*mm , // Target Thickness at 0 degree
- ThetaN );
-
- ThetaCM = myReaction -> EnergyLabToThetaCM( E , 1 ) /deg ;
- ResolThetaCM =ThetaCM - Init->GetICEmittedAngleThetaCM(0) ;
- Ex = myReaction -> ReconstructRelativistic( E , Theta ) ;
- X = HitDirection . X();
- Y = HitDirection . Y();
- }
+ E= He3TargetWind.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 15*micrometer , // Target Thickness at 0 degree
+ ThetaN );
+
+ E= He3TargetGaz.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 1.5*mm , // Target Thickness at 0 degree
+ ThetaN );
+
+ ThetaCM = myReaction -> EnergyLabToThetaCM( E , 1 ) /deg ;
+ ResolThetaCM =ThetaCM - Init->GetICEmittedAngleThetaCM(0) ;
+ Ex = myReaction -> ReconstructRelativistic( E , Theta ) ;
+ X = HitDirection . X();
+ Y = HitDirection . Y();
+ }
- else if(E>-1000 )
- {
- if(E>18)//CsI are inside a Mylar foil, plus rear alu strip
- {
- E= He3TargetWind.EvaluateInitialEnergy( E , // Energy of the detected particle
- 3*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
+ else if(E>-1000 )
+ {
+ if(E>18)//CsI are inside a Mylar foil, plus rear alu strip
+ {
+ // E= He3TargetWind.EvaluateInitialEnergy( E , // Energy of the detected particle
+ // 3*micrometer , // Target Thickness at 0 degree
+ // ThetaMM2Surface );
+ E= He3StripAl.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
+ }
+
E= He3StripAl.EvaluateInitialEnergy( E , // Energy of the detected particle
0.4*micrometer , // Target Thickness at 0 degree
ThetaMM2Surface );
+
+ E= He3TargetWind.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 15*micrometer , // Target Thickness at 0 degree
+ ThetaN );
+
+ E= He3TargetGaz.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 1.5*mm , // Target Thickness at 0 degree
+ ThetaN );
+ ThetaCM = myReaction -> EnergyLabToThetaCM( E , 1 ) /deg ;
+ Ex = myReaction -> ReconstructRelativistic( E, Theta ) ;
+ X = HitDirection . X();
+ Y = HitDirection . Y();
+
+ }
+
+ else {Ex=-100 ; X = -100 ; Y = -100 ;}
+
}
-
- E= He3StripAl.EvaluateInitialEnergy( E , // Energy of the detected particle
- 0.4*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
-
- E= He3TargetWind.EvaluateInitialEnergy( E , // Energy of the detected particle
- 15*micrometer , // Target Thickness at 0 degree
- ThetaN );
- E= He3TargetGaz.EvaluateInitialEnergy( E , // Energy of the detected particle
- 1.5*mm , // Target Thickness at 0 degree
- ThetaN );
+ if(M2 -> GetPositionOfInteraction().Z()<0)
+ {
+
+ if(E>-1000 )
+ {
+ if(E>18)
+ {
+ E= protonStripAl.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
+ }
+
+ E= protonStripAl.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
+
+ E= protonTargetWind.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 15*micrometer , // Target Thickness at 0 degree
+ ThetaN );
+
+ E= protonTargetGaz.EvaluateInitialEnergy( E , // Energy of the detected particle
+ 1.5*mm , // Target Thickness at 0 degree
+ ThetaN );
+ ThetaCM = myReaction -> EnergyLabToThetaCM( E , 1 ) /deg ;
+ Ex = myReaction -> ReconstructRelativistic( E, Theta ) ;
+ X = HitDirection . X();
+ Y = HitDirection . Y();
+
+ }
- Ex = myReaction -> ReconstructRelativistic( E, Theta ) ;
- X = HitDirection . X();
- Y = HitDirection . Y();
+ else {Ex=-100 ; X = -100 ; Y = -100 ;}
}
- else {Ex=-100 ; X = -100 ; Y = -100 ;}
EE = E ; TT = Theta/deg ;
diff --git a/NPLib/MUST2/Must2Array.cxx b/NPLib/MUST2/Must2Array.cxx
index ff4a1f8a3e6e9759c95f0cd2d15417922471f467..b627d4524f000f8bdea56abc9f951702fbae34bd 100644
--- a/NPLib/MUST2/Must2Array.cxx
+++ b/NPLib/MUST2/Must2Array.cxx
@@ -453,46 +453,31 @@ void MUST2Array::AddTelescope( double theta ,
TVector3 W ;
//Vector position of Telescope Face center
TVector3 C ;
-
- /*if(theta==180 && phi==90)
- {
- C = TVector3 (0,0,distance) ;
- U = TVector3 (1,0,0) ;
- V = TVector3 (0,1,0) ;
- W = TVector3 (0,0,1) ;
- }*/
-
- if(theta==0 && phi==0)
- {
- C = TVector3 (0,0,distance) ;
- U = TVector3 (1,0,0) ;
- V = TVector3 (0,1,0) ;
- W = TVector3 (0,0,1) ;
- }
- else
- {
- C = TVector3 ( distance * sin(theta) * cos(phi) ,
- distance * sin(theta) * sin(phi) ,
- distance * cos(theta) );
-
- W = C.Unit() ;
- U = W .Cross ( TVector3(0,1,0) ) ;
- V = W .Cross ( U );
-
- U = U.Unit();
- V = V.Unit();
-
- U.Rotate( beta_u * Pi/180. , U ) ;
- V.Rotate( beta_u * Pi/180. , U ) ;
-
- U.Rotate( beta_v * Pi/180. , V ) ;
- V.Rotate( beta_v * Pi/180. , V ) ;
-
- U.Rotate( beta_w * Pi/180. , W ) ;
- V.Rotate( beta_w * Pi/180. , W ) ;
- }
-
+ C = TVector3 ( distance * sin(theta) * cos(phi) ,
+ distance * sin(theta) * sin(phi) ,
+ distance * cos(theta) );
+
+ TVector3 P = TVector3( cos(theta ) * cos(phi) ,
+ cos(theta ) * sin(phi) ,
+ -sin(theta) );
+
+ W = C.Unit() ;
+ U = W .Cross ( P ) ;
+ V = W .Cross ( U );
+
+ U = U.Unit();
+ V = V.Unit();
+
+ U.Rotate( beta_u * Pi/180. , U ) ;
+ V.Rotate( beta_u * Pi/180. , U ) ;
+
+ U.Rotate( beta_v * Pi/180. , V ) ;
+ V.Rotate( beta_v * Pi/180. , V ) ;
+
+ U.Rotate( beta_w * Pi/180. , W ) ;
+ V.Rotate( beta_w * Pi/180. , W ) ;
+
double Face = 98 ; //mm
double NumberOfStrip = 128 ;
double StripPitch = Face/NumberOfStrip ; //mm
diff --git a/NPSimulation/src/GaspardTrackerSquare.cc b/NPSimulation/src/GaspardTrackerSquare.cc
index 16ccb0afd78b67a9a95378d2f27896133a255727..31f481245d2bf11fee67e971ec718337e29f12c4 100644
--- a/NPSimulation/src/GaspardTrackerSquare.cc
+++ b/NPSimulation/src/GaspardTrackerSquare.cc
@@ -760,11 +760,7 @@ void GaspardTrackerSquare::ConstructDetector(G4LogicalVolume* world)
else {
G4double Theta = m_Theta[i] ;
G4double Phi = m_Phi[i] ;
- //This part because if Phi and Theta = 0 equation are false
- if (Theta == 0) Theta = 0.0001 ;
- if (Theta == 2*cos(0)) Theta = 2 * acos(0) - 0.00001 ;
- if (Phi == 0) Phi = 0.0001 ;
-
+
// (u,v,w) unitary vector associated to telescope referencial
// (u,v) // to silicon plan
// w perpendicular to (u,v) plan and pointing ThirdStage
diff --git a/NPSimulation/src/MUST2Array.cc b/NPSimulation/src/MUST2Array.cc
index 1343f4ebb8b6b647aced63c522f3767f1f9fa22a..23cdf26108a9a66434b5e471a2ebe0139501df05 100644
--- a/NPSimulation/src/MUST2Array.cc
+++ b/NPSimulation/src/MUST2Array.cc
@@ -856,37 +856,32 @@ void MUST2Array::ConstructDetector(G4LogicalVolume* world)
else {
G4double Theta = m_Theta[i] ;
G4double Phi = m_Phi[i] ;
- //This part because if Phi and Theta = 0 equation are false
- if (Theta == 0) {
- Theta = 0.0001 ;
- }
- if (Theta == 2*acos(0)) {
- Theta = 2 * acos(0) - 0.00001 ;
- }
- if (Phi == 0) {
- Phi = 0.0001 ;
- }
-
+
// (u,v,w) unitary vector associated to telescope referencial
// (u,v) // to silicon plan
- // w perpendicular to (u,v) plan and pointing CsI
+ // w perpendicular to (u,v) plan and pointing ThirdStage
// Phi is angle between X axis and projection in (X,Y) plan
// Theta is angle between position vector and z axis
G4double wX = m_R[i] * sin(Theta / rad) * cos(Phi / rad) ;
G4double wY = m_R[i] * sin(Theta / rad) * sin(Phi / rad) ;
G4double wZ = m_R[i] * cos(Theta / rad) ;
-
MMw = G4ThreeVector(wX, wY, wZ) ;
- G4ThreeVector CT = MMw ;
- MMw = MMw.unit() ;
- G4ThreeVector Y = G4ThreeVector(0 , 1 , 0) ;
+ // vector corresponding to the center of the module
+ G4ThreeVector CT = MMw;
- MMu = MMw.cross(Y) ;
- MMv = MMw.cross(MMu) ;
+ // vector parallel to one axis of silicon plane
+ G4double ii = cos(Theta / rad) * cos(Phi / rad);
+ G4double jj = cos(Theta / rad) * sin(Phi / rad);
+ G4double kk = -sin(Theta / rad);
+ G4ThreeVector Y = G4ThreeVector(ii, jj, kk);
+ MMw = MMw.unit();
+ MMu = MMw.cross(Y);
+ MMv = MMw.cross(MMu);
MMv = MMv.unit();
MMu = MMu.unit();
+
// Passage Matrix from Lab Referential to Telescope Referential
// MUST2
MMrot = new G4RotationMatrix(MMu, MMv, MMw);
diff --git a/NPSimulation/src/ThinSi.cc b/NPSimulation/src/ThinSi.cc
index cd4eb8f94a1480765c1aee2b665e5d6f576094e8..494c4d01b027c517302211f00fb1a776c70a5144 100644
--- a/NPSimulation/src/ThinSi.cc
+++ b/NPSimulation/src/ThinSi.cc
@@ -558,24 +558,29 @@ void ThinSi::ConstructDetector(G4LogicalVolume* world)
// (u,v,w) unitary vector associated to telescope referencial
// (u,v) // to silicon plan
- // w perpendicular to (u,v) plan and pointing outside
+ // w perpendicular to (u,v) plan and pointing ThirdStage
// Phi is angle between X axis and projection in (X,Y) plan
// Theta is angle between position vector and z axis
G4double wX = m_R[i] * sin(Theta / rad) * cos(Phi / rad) ;
G4double wY = m_R[i] * sin(Theta / rad) * sin(Phi / rad) ;
G4double wZ = m_R[i] * cos(Theta / rad) ;
+ Det_w = G4ThreeVector(wX, wY, wZ) ;
- Det_w = G4ThreeVector(wX, wY, wZ) ;
- G4ThreeVector CT = Det_w ;
- Det_w = Det_w.unit() ;
-
- G4ThreeVector Y = G4ThreeVector(0 , 1 , 0) ;
+ // vector corresponding to the center of the module
+ G4ThreeVector CT = Det_w;
- Det_u = Det_w.cross(Y) ;
- Det_v = Det_w.cross(Det_u) ;
+ // vector parallel to one axis of silicon plane
+ G4double ii = cos(Theta / rad) * cos(Phi / rad);
+ G4double jj = cos(Theta / rad) * sin(Phi / rad);
+ G4double kk = -sin(Theta / rad);
+ G4ThreeVector Y = G4ThreeVector(ii, jj, kk);
+ Det_w = Det_w.unit();
+ Det_u = Det_w.cross(Y);
+ Det_v = Det_w.cross(Det_u);
Det_v = Det_v.unit();
Det_u = Det_u.unit();
+
// Passage Matrix from Lab Referential to Telescope Referential
// MUST2
Det_rot = new G4RotationMatrix(Det_u, Det_v, Det_w);
@@ -584,7 +589,7 @@ void ThinSi::ConstructDetector(G4LogicalVolume* world)
Det_rot->rotate(m_beta_v[i], Det_v);
Det_rot->rotate(m_beta_w[i], Det_w);
// translation to place Telescope
- Det_pos = CT ;
+ Det_pos = Det_w + CT ;
}