[G12] DOCA

[Michael C. Kunkel]

Greetings,

Here is the method for the DOCA


   vEg.SetPxPyPzE(0.0,0.0,egam,egam); //beam four-momemtum info
     vP.SetPxPyPzE(P_Px,P_Py,P_Pz,E_P);    //proton four-momemtum info
     vPip.SetPxPyPzE(Pip_Px,Pip_Py,Pip_Pz,E_Pip);        //pi+ 
four-momemtum info
     vPim.SetPxPyPzE(Pim_Px,Pim_Py,Pim_Pz,E_Pim);    //pi- four-momemtum 
info

//-----------------------------------------------DOCA-------------------------------------------------------------------------

     vertex3_P.SetXYZ(prot_vx, prot_vy, prot_vz);    //vertex of proton
     vertex3_Pip.SetXYZ(pip_vx, pip_vy, pip_vz);     //vertex of pi+
     vertex3_Pim.SetXYZ(pim_vx, pim_vy, pim_vz); //vertex of pi-
     vertex3_Eg.SetXYZ(prot_vx, prot_vy, prot_vz);  //vertex of beam 
photon this should be proton vertex for non-deteached vertices

     TVector3 PipPimMid = MatchTracks(vertex3_Pip, vPip.Vect(), 
vertex3_Pim, vPim.Vect(), &PipPimDist, &PipPimStat);

     PipPim_Midx = PipPimMid.X();
     PipPim_Midy = PipPimMid.Y();
     PipPim_Midz = PipPimMid.Z();

     TVector3 PPhotMid = MatchTracks(vertex3_P, vP.Vect(), vertex3_Eg, 
vEg.Vect(), &PPhotDist, &PPhotStat);

     TVector3 DIFFER =  PipPimMid - PPhotMid;
     DOCA_ang = (DIFFER.Dot(IV_PipPimVec.Vect()) / (DIFFER.Mag() * 
IV_PipPimVec.Vect().Mag()) ); //minimize this


The MatchTracks codes is as follows


TVector3  MatchTracks( TVector3 v1, TVector3 p1,  TVector3 v2, TVector3 
p2, Double_t *dist, Int_t *stat)
{
   Double_t R1 = 0;
   TVector3 diff = (v1) - (v2);

   TVector3 p1unit = p1.Unit();
   TVector3 p2unit = p2.Unit();

   Double_t R = p1unit.Dot(p2unit);

   if( R >= 1. ) {
     cout << "lines are parallel" << endl;
     (*stat) = -1;
     return(-1);
   }

   else {
     R1 = 1 / ( 1- R*R );
   }

   TVector3 RP21 = p2unit*R - p1unit;
   TVector3 RP12 = p2unit - p1unit*R;

   Double_t DOT21 = diff.Dot(RP21);
   Double_t DOT12 = diff.Dot(RP12);

   TVector3 M1 = v1 + p1unit*R1*DOT21;
   TVector3 M2 = v2 + p2unit*R1*DOT12;

   TVector3 diff1 = M2 - M1;
   if(M2.Mag() > M1.Mag()){
     (*dist) = -diff1.Mag();
   }
   else{(*dist) = diff1.Mag();}

   (*stat) = 1;
   TVector3 MidPoint = ( M1 + M2 ) * 0.5;
   return  MidPoint;
}


BR
MK
----------------------------------------
Michael C. Kunkel, PhD
Forschungszentrum Jülich
Nuclear Physics Institute and Juelich Center for Hadron Physics
Experimental Hadron Structure (IKP-1)
www.fz-juelich.de/ikp



-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.jlab.org/pipermail/g12/attachments/20160712/f74b0ac9/attachment.html>