diff --git a/NPLib/Detectors/Sofia/GladFieldMap.cxx b/NPLib/Detectors/Sofia/GladFieldMap.cxx index ac96df6dba96e202a8acbc852f290b58e01eedc3..ef72f6f80493f3af6d7593857d5b3c9af71faad3 100644 --- a/NPLib/Detectors/Sofia/GladFieldMap.cxx +++ b/NPLib/Detectors/Sofia/GladFieldMap.cxx @@ -92,7 +92,7 @@ double GladFieldMap::Delta(const double* parameter){ static vector<TVector3> pos; pos = Propagate(parameter[0], m_InitPos, m_InitDir, false); - pos.back().RotateY(-m_CentralTheta+m_Tilt); + //pos.back().RotateY(-m_CentralTheta+m_Tilt); diff = pos.back() - m_FinalPos; //cout << pos.back().X() << " " << pos.back().Z() << endl; @@ -110,10 +110,10 @@ double GladFieldMap::FindBrho(TVector3 Pos_init, TVector3 Dir_init, TVector3 Pos m_FinalPos = Pos_final; m_InitDir = m_InitDir.Unit(); - BrhoScan(6,11,1,m_InitPos, m_InitDir); + //BrhoScan(6,11,1,m_InitPos, m_InitDir); static double param[1]; param[0] = m_BrhoScan->Eval(m_FinalPos.X()); - return param[0]; + //return param[0]; m_min->Clear(); m_min->SetPrecision(1e-6); @@ -151,6 +151,28 @@ TGraph* GladFieldMap::BrhoScan(double Brho_min, double Brho_max, double Brho_ste return m_BrhoScan; } +////////////////////////////////////////////////////////////////////// +double GladFieldMap::GetFlightPath(TVector3 vStart, double Brho, TVector3 Pos, TVector3 Dir){ + + double FlightPath = 0; + + vector<TVector3> track; + track = Propagate(Brho, Pos, Dir, true); + + FlightPath += (Pos - vStart).Mag(); + + unsigned int vsize = track.size(); + for(unsigned int i=0; i<vsize-1; i++){ + TVector3 point1 = track[i]; + TVector3 point2 = track[i+1]; + + TVector3 v12 = point2 - point1; + FlightPath += v12.Mag(); + } + + return FlightPath; +} + ////////////////////////////////////////////////////////////////////// TVector3 GladFieldMap::PropagateToMWPC(TVector3 pos, TVector3 dir){ // go to MWPC3 frame reference diff --git a/NPLib/Detectors/Sofia/GladFieldMap.h b/NPLib/Detectors/Sofia/GladFieldMap.h index d9a956537e75ba270b0bbe714be286451810d2fe..375680d939b04f401903ce09b5d8f1447ce559bc 100644 --- a/NPLib/Detectors/Sofia/GladFieldMap.h +++ b/NPLib/Detectors/Sofia/GladFieldMap.h @@ -99,6 +99,7 @@ class GladFieldMap{ m_B = 2.2*m_Scale; } void SetBin(double val) {m_bin = val;} + void SetTimeStep(double val) {m_dt = val;} void SetCentralTheta(double val) {m_CentralTheta = val;} void Set_MWPC3_Position(double x, double y, double z) {m_MWPC3_POS = TVector3(x,y,z);} @@ -124,6 +125,7 @@ class GladFieldMap{ double GetZmax() {return m_z_max;} double GetCentralTheta() {return m_CentralTheta;} double GetBin() {return m_bin;} + double GetTimeStep() {return m_dt;} TVector3 Get_MWPC3_Position() {return m_MWPC3_POS;} public: @@ -137,6 +139,7 @@ class GladFieldMap{ TGraph* BrhoScan(double Brho_min, double Brho_max, double Brho_step, TVector3 pos, TVector3 dir); TVector3 CalculateIntersectionPoint(vector<TVector3> vPos); vector<TVector3> Propagate(double Brho, TVector3 Pos, TVector3 Dir, bool store); + double GetFlightPath(TVector3 vStart, double Brho, TVector3 Pos, TVector3 Dir); TVector3 PropagateToMWPC(TVector3 pos, TVector3 dir); void func(NPL::Particle& N, TVector3 Pos, TVector3 Imp, TVector3& xk, TVector3& pk); double FindBrho(TVector3 Pos_init, TVector3 Dir_init, TVector3 Pos_final); diff --git a/Projects/s455/Analysis.cxx b/Projects/s455/Analysis.cxx index b275e42e599d4715ec43f207a40cb989a8b9b12b..39b70b82f6348b4cb1a97c86eba13fd076a3d7e8 100644 --- a/Projects/s455/Analysis.cxx +++ b/Projects/s455/Analysis.cxx @@ -110,10 +110,12 @@ void Analysis::Init(){ m_GladField->SetGladTiltAngle(-14.*deg); m_GladField->LoadMap("GladFieldMap_50mm.dat"); m_GladField->SetBin(50); + m_GladField->SetTimeStep(0.8); m_GladField->SetCentralTheta(-20.*deg); //double Z_MWPC3 = 7.852*m; - double Z_MWPC3 = 3.402*m; + //double Z_MWPC3 = 3.402*m; + double Z_MWPC3 = 3.5*m; double X_MWPC3 = (Z_MWPC3 - m_GladField->GetGladTurningPoint().Z())*tan(m_GladField->GetCentralTheta()); m_GladField->Set_MWPC3_Position(X_MWPC3,0,Z_MWPC3); @@ -604,6 +606,7 @@ void Analysis::FissionFragmentAnalysis(){ double Z3lab = 0;; double Tilt = 14.*deg;//; TVector3 vZ = TVector3(0,0,1); + TVector3 vStart = TVector3(0,0,-DistanceStartToG); double XC = 0; double ZC = 0; double XB = 0; @@ -631,7 +634,7 @@ void Analysis::FissionFragmentAnalysis(){ X3lab = TofHit[i].x3*cos(Theta0) + XMW3; Z3lab = TofHit[i].x3*sin(Theta0) + ZMW3; - TVector3 vE = TVector3(X3lab, 0, Z3lab); + TVector3 vE = TVector3(X3lab, TofHit[i].y, Z3lab); TVector3 dir = TVector3(sin(TofHit[i].theta_in), 0, cos(TofHit[i].theta_in)); TofHit[i].x3lab = X3lab; TofHit[i].z3lab = Z3lab; @@ -670,13 +673,14 @@ void Analysis::FissionFragmentAnalysis(){ TVector3 v3lab = TVector3(X3lab,0,Z3lab); - TVector3 v1 = TVector3(XB,0,ZB); + TVector3 v1 = TVector3(XB,0,ZB)-vStart; TVector3 v3 = TVector3(X3lab-XD,0,Z3lab-ZD); TofHit[i].omega = abs(2.*asin(sqrt(pow(XD-XB,2) + pow(ZD-ZB,2))/(2*TofHit[i].rho))); double Path1 = v1.Mag(); double Path2 = TofHit[i].rho*TofHit[i].omega; double Path3 = v3.Mag(); - double PathLength = Path1 + Path2 + Path3 + 740. + 50.; + //double PathLength = Path1 + Path2 + Path3 + 740. + 50.; + double PathLength = m_GladField->GetFlightPath(vStart, TofHit[i].Brho, vA, dir) + 740. + 50; PathLength = PathLength/1000.; TofHit[i].flight_path = PathLength;