larutil::GeometryUtilities Class Reference

#include <GeometryUtilities.hh>

Inheritance diagram for larutil::GeometryUtilities:
larlight::larlight_base

List of all members.

Public Member Functions

void Reconfigure ()
Int_t Get3DaxisN (Int_t iplane0, Int_t iplane1, Double_t omega0, Double_t omega1, Double_t &phi, Double_t &theta) const
Double_t CalculatePitch (UInt_t iplane0, Double_t phi, Double_t theta) const
Double_t CalculatePitchPolar (UInt_t iplane0, Double_t phi, Double_t theta) const
Double_t Get3DSpecialCaseTheta (Int_t iplane0, Int_t iplane1, Double_t dw0, Double_t dw1) const
Double_t Get2Dangle (Double_t deltawire, Double_t deltatime) const
Double_t Get2Dangle (Double_t wireend, Double_t wirestart, Double_t timeend, Double_t timestart) const
double Get2Dangle (const larutil::PxPoint *endpoint, const larutil::PxPoint *startpoint) const
double Get2DangleFrom3D (unsigned int plane, double phi, double theta) const
double Get2DangleFrom3D (unsigned int plane, TVector3 dir_vector) const
Double_t Get2Dslope (Double_t deltawire, Double_t deltatime) const
Double_t Get2Dslope (Double_t wireend, Double_t wirestart, Double_t timeend, Double_t timestart) const
double Get2Dslope (const larutil::PxPoint *endpoint, const larutil::PxPoint *startpoint) const
Double_t Get2DDistance (Double_t wire1, Double_t time1, Double_t wire2, Double_t time2) const
double Get2DDistance (const larutil::PxPoint *point1, const larutil::PxPoint *point2) const
Double_t Get2DPitchDistance (Double_t angle, Double_t inwire, Double_t wire) const
Double_t Get2DPitchDistanceWSlope (Double_t slope, Double_t inwire, Double_t wire) const
Int_t GetPointOnLine (Double_t slope, Double_t intercept, Double_t wire1, Double_t time1, Double_t &wireout, Double_t &timeout) const
Int_t GetPointOnLine (Double_t slope, Double_t wirestart, Double_t timestart, Double_t wire1, Double_t time1, Double_t &wireout, Double_t &timeout) const
int GetPointOnLine (Double_t slope, const larutil::PxPoint *startpoint, const larutil::PxPoint *point1, larutil::PxPoint &pointout) const
int GetPointOnLine (double slope, double intercept, const larutil::PxPoint *point1, larutil::PxPoint &pointout) const
Int_t GetPointOnLineWSlopes (Double_t slope, Double_t intercept, Double_t ort_intercept, Double_t &wireout, Double_t &timeout) const
Int_t GetPointOnLineWSlopes (double slope, double intercept, double ort_intercept, larutil::PxPoint &pointonline) const
const larlight::hit * FindClosestHit (const std::vector< larlight::hit * > &hitlist, UInt_t wire, Double_t time) const
UInt_t FindClosestHitIndex (const std::vector< larlight::hit * > &hitlist, UInt_t wirein, Double_t timein) const
PxPoint Get2DPointProjection (Double_t *xyz, Int_t plane) const
PxPoint Get2DPointProjectionCM (std::vector< double > xyz, int plane) const
PxPoint Get2DPointProjectionCM (double *xyz, int plane) const
PxPoint Get2DPointProjectionCM (TLorentzVector *xyz, int plane) const
Double_t GetTimeTicks (Double_t x, Int_t plane) const
Int_t GetPlaneAndTPC (const larlight::hit *h, UInt_t &p, UInt_t &w) const
Int_t GetProjectedPoint (const PxPoint *p0, const PxPoint *p1, PxPoint &pN) const
Int_t GetYZ (const PxPoint *p0, const PxPoint *p1, Double_t *yz) const
Double_t PitchInView (UInt_t plane, Double_t phi, Double_t theta) const
void GetDirectionCosines (Double_t phi, Double_t theta, Double_t *dirs) const
void SelectLocalHitlist (const std::vector< larlight::hit * > &hitlist, std::vector< larlight::hit * > &hitlistlocal_index, Double_t wire_start, Double_t time_start, Double_t linearlimit, Double_t ortlimit, Double_t lineslopetest)
void SelectLocalHitlist (const std::vector< larlight::hit * > &hitlist, std::vector< UInt_t > &hitlistlocal_index, Double_t wire_start, Double_t time_start, Double_t linearlimit, Double_t ortlimit, Double_t lineslopetest)
void SelectLocalHitlist (const std::vector< larutil::PxHit > &hitlist, std::vector< const larutil::PxHit * > &hitlistlocal, larutil::PxPoint &startHit, Double_t &linearlimit, Double_t &ortlimit, Double_t &lineslopetest, larutil::PxHit &averageHit)
void SelectPolygonHitList (const std::vector< larutil::PxHit > &hitlist, std::vector< const larutil::PxHit * > &hitlistlocal)
std::vector< size_t > PolyOverlap (std::vector< const larutil::PxHit * > ordered_hits, std::vector< size_t > candidate_polygon)
bool Clockwise (double Ax, double Ay, double Bx, double By, double Cx, double Cy)
Double_t TimeToCm () const
Double_t WireToCm () const
Double_t WireTimeToCmCm () const
virtual void set_verbosity (MSG::Level level)
 Setter for the verbosity level.
MSG::Level get_verbosity () const
 Getter for the verbosity level.
const std::string class_name () const
 Getter for the class name.
void print (MSG::Level level, std::string where, std::string msg) const
 message print out method
void print (MSG::Level level, std::string msg) const
 message print out method

Static Public Member Functions

static const GeometryUtilities * GetME ()
 Singleton getter.

Protected Attributes

char _buf [200]
 char buffer for message manipulation
std::vector< bool > _verbosity
 holder for enabled message levels
MSG::Level _verbosity_level
 holder for specified verbosity level
std::string _name
 class name holder

Private Member Functions

 GeometryUtilities ()
 Default constructor = private for singleton.
 ~GeometryUtilities ()
 Default destructor.

Private Attributes

larutil::Geometry * geom
larutil::DetectorProperties * detp
larutil::LArProperties * larp
std::vector< Double_t > vertangle
Double_t fWirePitch
Double_t fTimeTick
Double_t fDriftVelocity
UInt_t fNPlanes
Double_t fWiretoCm
Double_t fTimetoCm
Double_t fWireTimetoCmCm

Static Private Attributes

static GeometryUtilities * _me = 0

Detailed Description

Definition at line 27 of file GeometryUtilities.hh.

Constructor & Destructor Documentation

larutil::GeometryUtilities::GeometryUtilities (  )  [private]

Default constructor = private for singleton.

Definition at line 17 of file GeometryUtilities.cc.

References larlight::larlight_base::_name, and Reconfigure().

00018   {
00019     _name = "GeometryUtilities";
00020 
00021     Reconfigure();
00022   }

larutil::GeometryUtilities::~GeometryUtilities (  )  [private]

Default destructor.

Definition at line 45 of file GeometryUtilities.cc.

00046   {
00047     
00048   }

Member Function Documentation

Double_t larutil::GeometryUtilities::CalculatePitch ( UInt_t  iplane0,
Double_t  phi,
Double_t  theta 
) const

Definition at line 294 of file GeometryUtilities.cc.

References larlight::MSG::ERROR, geom, larlight::GEO::k3D, larlight::GEO::kUnknown, larutil::Geometry::Nplanes(), larutil::Geometry::PlaneToView(), larlight::larlight_base::print(), larutil::Geometry::WireAngleToVertical(), and larutil::Geometry::WirePitch().

00297   {
00298     
00299     Double_t pitch = -1.;
00300     
00301     if(geom->PlaneToView(iplane) == larlight::GEO::kUnknown || 
00302        geom->PlaneToView(iplane) == larlight::GEO::k3D){
00303       print(larlight::MSG::ERROR,__FUNCTION__,
00304         Form("Warning :  no Pitch foreseen for view %d", geom->PlaneToView(iplane)));
00305       return pitch;
00306     }
00307     else{
00308      
00309       Double_t pi=TMath::Pi();
00310       
00311       Double_t fTheta=pi/2-theta;
00312      
00313       Double_t fPhi=-(phi+pi/2);
00314       //Double_t fPhi=pi/2-phi;
00315       //if(fPhi<0)
00316       //    fPhi=phi-pi/2;
00317      
00318       //fTheta=TMath::Pi()/2;
00319  
00320      
00321      
00322       for(UInt_t i = 0; i < geom->Nplanes(); ++i) {
00323     if(i == iplane){
00324       Double_t wirePitch = geom->WirePitch(0,1,i);
00325       Double_t angleToVert =0.5*TMath::Pi() - geom->WireAngleToVertical(geom->PlaneToView(i));
00326       
00327           //    //    //std::cout <<" %%%%%%%%%%  " << i << " angle " 
00328           //                       << angleToVert*180/pi << " " 
00329           //                       << geom->Plane(i).Wire(0).ThetaZ(false)*180/pi 
00330           //                       <<" wirePitch " << wirePitch
00331           //                       <<"\n %%%%%%%%%%  " << fTheta << " " << fPhi<< std::endl;
00332           //       
00333       
00334       Double_t cosgamma = TMath::Abs(TMath::Sin(angleToVert)*TMath::Cos(fTheta)
00335                        +TMath::Cos(angleToVert)*TMath::Sin(fTheta)*TMath::Sin(fPhi));
00336       
00337       if (cosgamma>0) pitch = wirePitch/cosgamma;     
00338     } // end if the correct view
00339       } // end loop over planes
00340     } // end if a reasonable view
00341    
00342     return pitch;
00343   }

Double_t larutil::GeometryUtilities::CalculatePitchPolar ( UInt_t  iplane0,
Double_t  phi,
Double_t  theta 
) const

Definition at line 351 of file GeometryUtilities.cc.

References larlight::MSG::ERROR, geom, larlight::GEO::k3D, larlight::GEO::kUnknown, larutil::Geometry::Nplanes(), larutil::Geometry::PlaneToView(), larlight::larlight_base::print(), larutil::Geometry::WireAngleToVertical(), and larutil::Geometry::WirePitch().

00354   {
00355 
00356     Double_t pitch = -1.;
00357   
00358     if(geom->PlaneToView(iplane) == larlight::GEO::kUnknown || 
00359        geom->PlaneToView(iplane) == larlight::GEO::k3D){
00360       print(larlight::MSG::ERROR,__FUNCTION__,
00361         Form("Warning :  no Pitch foreseen for view %d", geom->PlaneToView(iplane)));
00362       return pitch;
00363     }
00364     else{
00365         
00366       Double_t fTheta=theta;
00367       Double_t fPhi=phi;  
00368      
00369      
00370       //fTheta=TMath::Pi()/2;
00371      
00372      
00373      
00374       for(UInt_t i = 0; i < geom->Nplanes(); ++i){
00375     if(i == iplane){
00376       Double_t wirePitch = geom->WirePitch(0,1,i);
00377       Double_t angleToVert =0.5*TMath::Pi() - geom->WireAngleToVertical(geom->PlaneToView(i));
00378       
00379       //        //std::cout <<" %%%%%%%%%%  " << i << " angle " 
00380       //                       << angleToVert*180/pi << " " 
00381       //                       << geom->Plane(i).Wire(0).ThetaZ(false)*180/pi 
00382       //                       <<" wirePitch " << wirePitch
00383       //                       <<"\n %%%%%%%%%%  " << fTheta << " " << fPhi<< std::endl;
00384       
00385       
00386       Double_t cosgamma = TMath::Abs(TMath::Sin(angleToVert)*TMath::Cos(fTheta)
00387                        +TMath::Cos(angleToVert)*TMath::Sin(fTheta)*TMath::Sin(fPhi));
00388       
00389       if (cosgamma>0) pitch = wirePitch/cosgamma;     
00390     } // end if the correct view
00391       } // end loop over planes
00392     } // end if a reasonable view
00393    
00394     return pitch;
00395   }

const std::string larlight::larlight_base::class_name (  )  const [inline, inherited]

Getter for the class name.

Definition at line 49 of file larlight_base.hh.

References larlight::larlight_base::_name.

00049 {return _name;};

bool larutil::GeometryUtilities::Clockwise ( double  Ax,
double  Ay,
double  Bx,
double  By,
double  Cx,
double  Cy 
)

Definition at line 1348 of file GeometryUtilities.cc.

Referenced by PolyOverlap().

01348                                                                                               {
01349     return (Cy-Ay)*(Bx-Ax) > (By-Ay)*(Cx-Ax);
01350   }

const larlight::hit * larutil::GeometryUtilities::FindClosestHit ( const std::vector< larlight::hit * > &  hitlist,
UInt_t  wire,
Double_t  time 
) const

Definition at line 742 of file GeometryUtilities.cc.

References FindClosestHitIndex().

00745   {
00746     return hitlist.at(FindClosestHitIndex(hitlist,wirein,timein));
00747   }

UInt_t larutil::GeometryUtilities::FindClosestHitIndex ( const std::vector< larlight::hit * > &  hitlist,
UInt_t  wirein,
Double_t  timein 
) const

Definition at line 753 of file GeometryUtilities.cc.

References Get2DDistance(), GetPlaneAndTPC(), larlight::DATA::INVALID_DOUBLE, and larlight::DATA::INVALID_UINT.

Referenced by FindClosestHit().

00756   {
00757   
00758     Double_t min_length_from_start=larlight::DATA::INVALID_DOUBLE;
00759     UInt_t index=larlight::DATA::INVALID_UINT;
00760    
00761     UInt_t plane,wire;
00762    
00763    
00764     for(UInt_t ii=0; ii<hitlist.size();ii++){
00765 
00766       Double_t time = hitlist.at(ii)->PeakTime();  
00767       GetPlaneAndTPC(hitlist.at(ii),plane,wire);
00768       
00769       Double_t dist_mod=Get2DDistance(wirein,timein,wire,time);
00770 
00771       if(dist_mod<min_length_from_start){
00772     //wire_start[plane]=wire;
00773     //time_start[plane]=time;
00774     index = ii;
00775     min_length_from_start=dist_mod;
00776       } 
00777       
00778     } 
00779   
00780     return index;
00781   }

Double_t larutil::GeometryUtilities::Get2Dangle ( const larutil::PxPoint *  endpoint,
const larutil::PxPoint *  startpoint 
) const

Definition at line 459 of file GeometryUtilities.cc.

References Get2Dangle().

00461   {
00462     return Get2Dangle(endpoint->w - startpoint->w, endpoint->t - startpoint->t);
00463   
00464   }

Double_t larutil::GeometryUtilities::Get2Dangle ( Double_t  wireend,
Double_t  wirestart,
Double_t  timeend,
Double_t  timestart 
) const

Definition at line 445 of file GeometryUtilities.cc.

References fTimetoCm, fWiretoCm, and Get2Dangle().

00449   {
00450 
00451     return Get2Dangle((wireend-wirestart)*fWiretoCm,(timeend-timestart)*fTimetoCm);
00452   
00453   }

Double_t larutil::GeometryUtilities::Get2Dangle ( Double_t  deltawire,
Double_t  deltatime 
) const

Definition at line 470 of file GeometryUtilities.cc.

Referenced by Get2Dangle(), Get2DangleFrom3D(), and Get2Dslope().

00472   {
00473  
00474       
00475     Double_t BC,AC;
00476     Double_t omega;
00477  
00478     BC = ((Double_t)dwire); // in cm
00479     AC = ((Double_t)dtime); //in cm 
00480     omega = std::asin(  AC/std::sqrt(pow(AC,2)+pow(BC,2)) );
00481     if(BC<0)  // for the time being. Will check if it works for AC<0
00482       { 
00483     if(AC>0){
00484       omega= TMath::Pi()-std::abs(omega);  //
00485     }
00486     else if(AC<0){
00487       omega=-TMath::Pi()+std::abs(omega);
00488     }
00489     else {
00490       omega=TMath::Pi();
00491     }
00492       } 
00493     //return omega;
00494     return omega; //*fWirePitch/(fTimeTick*fDriftVelocity);
00495 
00496   }

double larutil::GeometryUtilities::Get2DangleFrom3D ( unsigned int  plane,
TVector3  dir_vector 
) const

Definition at line 513 of file GeometryUtilities.cc.

References larutil::Geometry::DetHalfHeight(), larutil::Geometry::DetHalfWidth(), larutil::Geometry::DetLength(), geom, Get2Dangle(), larutil::Geometry::PlaneToView(), and larutil::Geometry::WireAngleToVertical().

00514   {
00515    double alpha= 0.5*TMath::Pi()-geom->WireAngleToVertical(geom->PlaneToView(plane)); 
00516    // create dummy  xyz point in middle of detector and another one in unit length.
00517    // calculate correspoding points in wire-time space and use the differnces between those to return 2D a
00518    // angle
00519    
00520    
00521    TVector3 start(geom->DetHalfWidth(),0.,geom->DetLength()/2.);
00522    TVector3 end=start+dir_vector;
00523    
00524    
00525     //the wire coordinate is already in cm. The time needs to be converted.
00526    larutil::PxPoint startp(plane,(geom->DetHalfHeight()*sin(fabs(alpha))+start[2]*cos(alpha)-start[1]*sin(alpha)),start[0]);
00527    
00528    larutil::PxPoint endp(plane,(geom->DetHalfHeight()*sin(fabs(alpha))+end[2]*cos(alpha)-end[1]*sin(alpha)),end[0]);
00529    
00530    double angle=Get2Dangle(&endp,&startp);
00531    
00532       
00533    return angle;
00534     
00535   }

double larutil::GeometryUtilities::Get2DangleFrom3D ( unsigned int  plane,
double  phi,
double  theta 
) const

Definition at line 501 of file GeometryUtilities.cc.

00502    {
00503    TVector3 dummyvector(cos(theta*TMath::Pi()/180.)*sin(phi*TMath::Pi()/180.)  ,sin(theta*TMath::Pi()/180.)  ,  cos(theta*TMath::Pi()/180.)*cos(phi*TMath::Pi()/180.));
00504     
00505     return Get2DangleFrom3D(plane,dummyvector);
00506      
00507    }

double larutil::GeometryUtilities::Get2DDistance ( const larutil::PxPoint *  point1,
const larutil::PxPoint *  point2 
) const

Definition at line 556 of file GeometryUtilities.cc.

00558   {
00559     return TMath::Sqrt( pow((point1->w - point2->w),2)+pow((point1->t - point2->t),2) );    
00560   }

Double_t larutil::GeometryUtilities::Get2DDistance ( Double_t  wire1,
Double_t  time1,
Double_t  wire2,
Double_t  time2 
) const

Definition at line 545 of file GeometryUtilities.cc.

References fTimetoCm, and fWiretoCm.

Referenced by FindClosestHitIndex(), and SelectLocalHitlist().

00549   {
00550     
00551     return TMath::Sqrt( pow((wire1-wire2)*fWiretoCm,2)+pow((time1-time2)*fTimetoCm,2) );    
00552   
00553   }

Double_t larutil::GeometryUtilities::Get2DPitchDistance ( Double_t  angle,
Double_t  inwire,
Double_t  wire 
) const

Definition at line 567 of file GeometryUtilities.cc.

References fWiretoCm.

00570   {
00571     Double_t radangle = TMath::Pi()*angle/180;
00572     if(std::cos(radangle)==0)
00573       return 9999;
00574     return std::abs((wire-inwire)/std::cos(radangle))*fWiretoCm; 
00575   }

Double_t larutil::GeometryUtilities::Get2DPitchDistanceWSlope ( Double_t  slope,
Double_t  inwire,
Double_t  wire 
) const

Definition at line 579 of file GeometryUtilities.cc.

References fWiretoCm.

00582   {
00583   
00584     return std::abs(wire-inwire)*TMath::Sqrt(1+slope*slope)*fWiretoCm; 
00585    
00586   }

PxPoint larutil::GeometryUtilities::Get2DPointProjection ( Double_t *  xyz,
Int_t  plane 
) const

Todo:
: this should use the cryostat and tpc as well in the NearestWire method

Definition at line 881 of file GeometryUtilities.cc.

References detp, fDriftVelocity, fTimeTick, geom, larutil::Geometry::NearestWire(), larutil::Geometry::PlaneOriginVtx(), and larutil::DetectorProperties::TriggerOffset().

Referenced by GetProjectedPoint().

00881                                                                                  {
00882   
00883     PxPoint pN(0,0,0);
00884     
00885     Double_t pos[3];
00886     geom->PlaneOriginVtx(plane,pos);
00887     Double_t drifttick=(xyz[0]/fDriftVelocity)*(1./fTimeTick);
00888       
00889     pos[1]=xyz[1];
00890     pos[2]=xyz[2];
00891 
00893     
00894     pN.w = geom->NearestWire(pos, plane);
00895     pN.t=drifttick-(pos[0]/fDriftVelocity)*(1./fTimeTick)+detp->TriggerOffset();  
00896     pN.plane=plane;
00897     
00898     return pN;
00899      
00900    }

PxPoint larutil::GeometryUtilities::Get2DPointProjectionCM ( TLorentzVector *  xyz,
int  plane 
) const

Definition at line 952 of file GeometryUtilities.cc.

References Get2DPointProjectionCM().

00952                                                                                        {
00953    double xyznew[3]={(*xyz)[0],(*xyz)[1],(*xyz)[2]}; 
00954     
00955    return  Get2DPointProjectionCM(xyznew,plane);
00956   }

PxPoint larutil::GeometryUtilities::Get2DPointProjectionCM ( double *  xyz,
int  plane 
) const

Todo:
: this should use the cryostat and tpc as well in the NearestWire method

Definition at line 930 of file GeometryUtilities.cc.

References fWiretoCm, geom, and larutil::Geometry::NearestWire().

00930                                                                                {
00931   
00932     PxPoint pN(0,0,0);
00933     
00934     Double_t pos[3];
00935         
00936       
00937     pos[1]=xyz[1];
00938     pos[2]=xyz[2];
00939 
00941     
00942     pN.w = geom->NearestWire(pos, plane)*fWiretoCm;
00943     pN.t=xyz[0];  
00944     pN.plane=plane;
00945     
00946     return pN;
00947      
00948    }

PxPoint larutil::GeometryUtilities::Get2DPointProjectionCM ( std::vector< double >  xyz,
int  plane 
) const

Todo:
: this should use the cryostat and tpc as well in the NearestWire method

Definition at line 909 of file GeometryUtilities.cc.

References fWiretoCm, geom, and larutil::Geometry::NearestWire().

Referenced by Get2DPointProjectionCM().

00909                                                                                            {
00910   
00911     PxPoint pN(0,0,0);
00912     
00913     Double_t pos[3];
00914         
00915       
00916     pos[1]=xyz[1];
00917     pos[2]=xyz[2];
00918 
00920     
00921     pN.w = geom->NearestWire(pos, plane)*fWiretoCm;
00922     pN.t=xyz[0];  
00923     pN.plane=plane;
00924     
00925     return pN;
00926      
00927    }

double larutil::GeometryUtilities::Get2Dslope ( const larutil::PxPoint *  endpoint,
const larutil::PxPoint *  startpoint 
) const

Definition at line 418 of file GeometryUtilities.cc.

References Get2Dslope().

00420   {
00421     return Get2Dslope(endpoint->w - startpoint->w,endpoint->t - startpoint->t);
00422   
00423   }

Double_t larutil::GeometryUtilities::Get2Dslope ( Double_t  wireend,
Double_t  wirestart,
Double_t  timeend,
Double_t  timestart 
) const

Definition at line 404 of file GeometryUtilities.cc.

References fTimetoCm, fWiretoCm, and Get2Dslope().

00408   {
00409     
00410     return GeometryUtilities::Get2Dslope((wireend-wirestart)*fWiretoCm,(timeend-timestart)*fTimetoCm);
00411   
00412   }

Double_t larutil::GeometryUtilities::Get2Dslope ( Double_t  deltawire,
Double_t  deltatime 
) const

Definition at line 430 of file GeometryUtilities.cc.

References fWireTimetoCmCm, and Get2Dangle().

Referenced by Get2Dslope().

00432   {
00433  
00434     //return omega;
00435  
00436     return tan(Get2Dangle(dwire,dtime))/fWireTimetoCmCm;
00437 
00438   }

Int_t larutil::GeometryUtilities::Get3DaxisN ( Int_t  iplane0,
Int_t  iplane1,
Double_t  omega0,
Double_t  omega1,
Double_t &  phi,
Double_t &  theta 
) const

Definition at line 58 of file GeometryUtilities.cc.

References larlight::DATA::INVALID_DOUBLE, and vertangle.

00064   {
00065  
00066     Double_t l(0),m(0),n(0);
00067     Double_t ln(0),mn(0),nn(0);
00068     Double_t phis(0),thetan(0);
00069     //Double_t phin(0);//,phis(0),thetan(0);
00070     // pretend collection and induction planes. 
00071     // "Collection" is the plane with the vertical angle equal to zero. 
00072     // If both are non zero collection is the one with the negative angle. 
00073     UInt_t Cplane=0,Iplane=1;   
00074     //then angleC and angleI are the respective angles to vertical in these planes and slopeC, slopeI are the tangents of the track.
00075     Double_t angleC,angleI,slopeC,slopeI,omegaC,omegaI;
00076     omegaC = larlight::DATA::INVALID_DOUBLE;
00077     omegaI = larlight::DATA::INVALID_DOUBLE;
00078     // don't know how to reconstruct these yet, so exit with error.
00079     
00080     if(omega0==0 || omega1==0){
00081       phi=0;
00082       theta=-999;
00083       return -1;
00084     }
00085     
00086     
00088     
00089     //check if backwards going track
00090     Double_t backwards=0;
00091     
00092     Double_t alt_backwards=0;
00093     
00095     if(fabs(omega0)>(TMath::Pi()/2.0) && fabs(omega1)>(TMath::Pi()/2.0) ) {
00096       backwards=1;
00097     }
00098     
00099     if(fabs(omega0)>(TMath::Pi()/2.0) || fabs(omega1)>(TMath::Pi()/2.0) ) {
00100       alt_backwards=1;
00101     }
00102     
00103     
00104     
00105     
00106     if(vertangle[iplane0] == 0){   
00107       // first plane is at 0 degrees
00108       Cplane=iplane0;
00109       Iplane=iplane1;
00110       omegaC=omega0;
00111       omegaI=omega1;
00112     }
00113     else if(vertangle[iplane1] == 0){  
00114       // second plane is at 0 degrees
00115       Cplane = iplane1;
00116       Iplane = iplane0;
00117       omegaC=omega1;
00118       omegaI=omega0;
00119     }
00120     else if(vertangle[iplane0] != 0 && vertangle[iplane1] != 0){
00121       //both planes are at non zero degree - find the one with deg<0
00122       if(vertangle[iplane1] < vertangle[iplane0]){
00123     Cplane = iplane1;
00124     Iplane = iplane0;
00125     omegaC=omega1;
00126     omegaI=omega0;
00127       }
00128       else if(vertangle[iplane1] > vertangle[iplane0]){
00129     Cplane = iplane0;
00130     Iplane = iplane1;
00131     omegaC=omega0;
00132     omegaI=omega1;
00133       }
00134       else{
00135     //throw error - same plane.
00136     return -1;
00137       } 
00138       
00139     }
00140     slopeC=tan(omegaC);
00141     slopeI=tan(omegaI);
00142     //omega0=tan(omega0);
00143     //omega1=tan(omega1);
00144     angleC=vertangle[Cplane];
00145     angleI=vertangle[Iplane];
00146     
00147     //0 -1 factor depending on if one of the planes is vertical.
00148     bool nfact = !(vertangle[Cplane]);
00149     
00150     
00151     
00152     if(omegaC < TMath::Pi() && omegaC > 0 )
00153       ln=1;
00154     else
00155       ln=-1;
00156     
00157     //std::cout << " slopes, C:"<< slopeC << " " << (omegaC) << " I:" << slopeI << " " << omegaI <<std::endl;
00158     slopeI=tan(omegaI);
00159     
00160     //std::cout << "omegaC, angleC " << omegaC << " " << angleC << "cond: " << omegaC-angleC << " ln: " << ln << std::endl;
00161     
00162     l = 1;
00163     
00164     
00165     m = (1/(2*sin(angleI)))*((cos(angleI)/(slopeC*cos(angleC)))-(1/slopeI) 
00166                  +nfact*(  cos(angleI)/slopeC-1/slopeI  )     );
00167     
00168     n = (1/(2*cos(angleC)))*((1/slopeC)+(1/slopeI) +nfact*((1/slopeC)-(1/slopeI)));
00169     
00170     mn = (ln/(2*sin(angleI)))*((cos(angleI)/(slopeC*cos(angleC)))-(1/slopeI) 
00171                    +nfact*(  cos(angleI)/(cos(angleC)*slopeC)-1/slopeI  )     );
00172     
00173     nn = (ln/(2*cos(angleC)))*((1/slopeC)+(1/slopeI) +nfact*((1/slopeC)-(1/slopeI)));
00174     
00175     
00176     float Phi;
00177     float alt_Phi;
00178     // Direction angles
00179     if(fabs(angleC)>0.01)  // catch numeric error values 
00180       {
00181     phi=atan(n/l);
00182     //phin=atan(ln/nn);
00183     phis=asin(ln/TMath::Sqrt(ln*ln+nn*nn));
00184         
00185     if(fabs(slopeC+slopeI) < 0.001)
00186       phis=0;
00187     else if(fabs(omegaI)>0.01 && (omegaI/fabs(omegaI) == -omegaC/fabs(omegaC) ) && ( fabs(omegaC) < 20*TMath::Pi()/180 || fabs(omegaC) > 160*TMath::Pi()/180   ) ) // angles have 
00188       {phis = (fabs(omegaC) > TMath::Pi()/2) ? TMath::Pi() : 0;    //angles are 
00189         
00190       }
00191     
00192     
00193     
00194     if(nn<0 && phis>0 && !(!alt_backwards && fabs(phis)<TMath::Pi()/4 ) )   // do not go back if track looks forward and phi is forward
00195       phis=(TMath::Pi())-phis;
00196     else if(nn<0 && phis<0 && !(!alt_backwards && fabs(phis)<TMath::Pi()/4 ) )
00197       phis=(-TMath::Pi())-phis;
00198     
00199     
00200     // solve the ambiguities due to tangent periodicity
00201     Phi = phi > 0. ? (TMath::Pi()/2)-phi : fabs(phi)-(TMath::Pi()/2) ; 
00202     alt_Phi = phi > 0. ? (TMath::Pi()/2)-phi : fabs(phi)-(TMath::Pi()/2) ; 
00203     
00204     if(backwards==1){
00205       if(Phi<0){ Phi=Phi+TMath::Pi();}
00206       else if(Phi>0){Phi=Phi-TMath::Pi();}
00207     }
00208     
00209     bool alt_condition=( ( fabs(omegaC)>0.75*TMath::Pi() && fabs(omegaI)>0.166*TMath::Pi() )|| ( fabs(omegaI)>0.75*TMath::Pi() && fabs(omegaC)>0.166*TMath::Pi() ) );
00210     
00211     
00212     if((alt_backwards==1 && alt_condition)   || backwards==1 ){
00213       if(alt_Phi<0){alt_Phi=alt_Phi+TMath::Pi();}
00214       else if(alt_Phi>0){alt_Phi=alt_Phi-TMath::Pi();}
00215     }
00216     
00217       }
00218     else  // if plane is collection than Phi = omega
00219       {phi=omegaC;
00220     Phi=omegaC;
00221     phis=omegaC;
00222     alt_Phi=omegaC;
00223       }
00224     
00225     
00226     theta = acos( m / (sqrt(pow(l,2)+pow(m,2)+pow(n,2)) ) ) ;
00227     thetan = -asin ( mn / (sqrt(pow(l,2)+pow(mn,2)+pow(nn,2)) ) ) ;
00228     //Double_t thetah = acos( mn / (sqrt(pow(l,2)+pow(mn,2)+pow(nn,2)) ) ) ;
00229     //float Theta;
00230     //float alt_Theta = 0.;
00231     
00232     
00233     
00234     
00235     //if(Phi < 0)Theta = (TMath::Pi()/2)-theta;
00236     //if(Phi > 0)Theta = theta-(TMath::Pi()/2);
00237     
00238     //if(alt_Phi < 0)alt_Theta = (TMath::Pi()/2)-theta;
00239     //if(alt_Phi > 0)alt_Theta = theta-(TMath::Pi()/2);
00240     
00242     //std::cout << "++++++++ GeomUtil_angles: Phi: " << alt_Phi*180/TMath::Pi() << " Theta: " << alt_Theta*180/TMath::Pi() << std::endl;
00243     
00244     //std::cout << "++++++++ GeomUtil_new_angles: Phi: " << phis*180/TMath::Pi() << " Theta: " << thetan*180/TMath::Pi() << std::endl;
00245     
00246     phi=phis*180/TMath::Pi();
00247     theta=thetan*180/TMath::Pi();
00248     
00249     
00250     return 0;   }

Double_t larutil::GeometryUtilities::Get3DSpecialCaseTheta ( Int_t  iplane0,
Int_t  iplane1,
Double_t  dw0,
Double_t  dw1 
) const

Definition at line 256 of file GeometryUtilities.cc.

References vertangle.

00260   {
00261   
00262     Double_t splane,lplane;   // plane in which the track is shorter and longer.
00263     Double_t sdw,ldw; 
00264   
00265     if(dw0==0 && dw1==0)
00266       return -1;
00267   
00268     if(dw0 >= dw1 ) {
00269       lplane=iplane0; 
00270       splane=iplane1;
00271       ldw=dw0;
00272       sdw=dw1;
00273     }
00274     else {
00275       lplane=iplane1; 
00276       splane=iplane0;
00277       ldw=dw1;
00278       sdw=dw0;
00279     }
00280   
00281     Double_t top=(std::cos(vertangle[splane])-std::cos(vertangle[lplane])*sdw/ldw);
00282     Double_t bottom = tan(vertangle[lplane]*std::cos(vertangle[splane]) ); 
00283           bottom -= tan(vertangle[splane]*std::cos(vertangle[lplane]) )*sdw/ldw;
00284   
00285     Double_t tantheta=top/bottom;
00286   
00287     return atan(tantheta)*vertangle[lplane]/std::abs(vertangle[lplane])*180./TMath::Pi();
00288   }

MSG::Level larlight::larlight_base::get_verbosity (  )  const [inline, inherited]

Getter for the verbosity level.

Definition at line 46 of file larlight_base.hh.

References larlight::larlight_base::_verbosity_level.

00046 {return _verbosity_level;};

void larutil::GeometryUtilities::GetDirectionCosines ( Double_t  phi,
Double_t  theta,
Double_t *  dirs 
) const

Definition at line 1008 of file GeometryUtilities.cc.

Referenced by PitchInView().

01011   {
01012     theta*=(TMath::Pi()/180);
01013     phi*=(TMath::Pi()/180); // working on copies, it's ok.
01014     dirs[0]=TMath::Cos(theta)*TMath::Sin(phi);
01015     dirs[1]=TMath::Sin(theta);
01016     dirs[2]=TMath::Cos(theta)*TMath::Cos(phi);
01017    
01018   }

static const GeometryUtilities* larutil::GeometryUtilities::GetME (  )  [inline, static]

Singleton getter.

Definition at line 41 of file GeometryUtilities.hh.

References _me.

Referenced by larutil::LArUtilManager::ReconfigureUtilities().

00041                                            {
00042       if(!_me) _me = new GeometryUtilities;
00043       return _me;
00044     }

Int_t larutil::GeometryUtilities::GetPlaneAndTPC ( const larlight::hit *  h,
UInt_t &  p,
UInt_t &  w 
) const

Definition at line 1022 of file GeometryUtilities.cc.

References larlight::hit::Channel(), larutil::Geometry::ChannelToPlane(), larutil::Geometry::ChannelToWire(), and geom.

Referenced by FindClosestHitIndex(), and SelectLocalHitlist().

01025   {
01026     p  = geom->ChannelToPlane(h->Channel());
01027     w  = geom->ChannelToWire(h->Channel());
01028     return 0;
01029   }

int larutil::GeometryUtilities::GetPointOnLine ( double  slope,
double  intercept,
const larutil::PxPoint *  point1,
larutil::PxPoint &  pointout 
) const

Definition at line 642 of file GeometryUtilities.cc.

00646   {
00647     double invslope=0;
00648       
00649     if(slope)   
00650     {
00651 //  invslope=-1./slope*fWireTimetoCmCm*fWireTimetoCmCm;
00652        invslope=-1./slope;
00653     }
00654   
00655     double ort_intercept=point1->t - invslope * point1->w;
00656     
00657     if((slope-invslope)!=0)
00658       pointout.w = (ort_intercept - intercept)/(slope-invslope); 
00659     else
00660       pointout.w = point1->w;
00661     
00662     pointout.t = slope * pointout.w + intercept;   
00663     
00664     return 0;
00665   }

int larutil::GeometryUtilities::GetPointOnLine ( Double_t  slope,
const larutil::PxPoint *  startpoint,
const larutil::PxPoint *  point1,
larutil::PxPoint &  pointout 
) const

Definition at line 624 of file GeometryUtilities.cc.

References GetPointOnLine().

00628   {
00629     
00630     double intercept=startpoint->t - slope * startpoint->w;  
00631     
00632     
00633     return GetPointOnLine(slope,intercept,point1,pointout);
00634     
00635   }

Int_t larutil::GeometryUtilities::GetPointOnLine ( Double_t  slope,
Double_t  wirestart,
Double_t  timestart,
Double_t  wire1,
Double_t  time1,
Double_t &  wireout,
Double_t &  timeout 
) const

Definition at line 671 of file GeometryUtilities.cc.

References GetPointOnLine().

00678   {
00679     Double_t intercept=timestart-slope*(Double_t)wirestart;
00680   
00681     return GetPointOnLine(slope,intercept,wire1,time1,wireout,timeout);
00682   }

Int_t larutil::GeometryUtilities::GetPointOnLine ( Double_t  slope,
Double_t  intercept,
Double_t  wire1,
Double_t  time1,
Double_t &  wireout,
Double_t &  timeout 
) const

Definition at line 595 of file GeometryUtilities.cc.

Referenced by GetPointOnLine(), and SelectLocalHitlist().

00601   {
00602     Double_t invslope=0;
00603       
00604     if(slope)   
00605       {
00606     invslope=-1./slope;
00607       }
00608   
00609     Double_t ort_intercept=time1-invslope*(Double_t)wire1;
00610     
00611     if((slope-invslope)!=0)
00612       wireout=(ort_intercept - intercept)/(slope-invslope); 
00613     else
00614       wireout=wire1;
00615     timeout=slope*wireout+intercept;   
00616     
00617     return 0;
00618   }

Int_t larutil::GeometryUtilities::GetPointOnLineWSlopes ( double  slope,
double  intercept,
double  ort_intercept,
larutil::PxPoint &  pointonline 
) const

Definition at line 718 of file GeometryUtilities.cc.

00722   {
00723     Double_t invslope=0;
00724   
00725     if(slope)   
00726       invslope=-1./slope;
00727 
00728 
00729     // invslope *= fWireTimetoCmCm * fWireTimetoCmCm;
00730   
00731     pointonline.w = (ort_intercept - intercept)/(slope-invslope); 
00732     pointonline.t = slope * pointonline.w + intercept; 
00733     return 0;  
00734   }

Int_t larutil::GeometryUtilities::GetPointOnLineWSlopes ( Double_t  slope,
Double_t  intercept,
Double_t  ort_intercept,
Double_t &  wireout,
Double_t &  timeout 
) const

Definition at line 688 of file GeometryUtilities.cc.

References fTimetoCm, fWireTimetoCmCm, and fWiretoCm.

00693   {
00694     Double_t invslope=0;
00695   
00696     if(slope)   
00697     {
00698         invslope=-1./slope;
00699     }
00700     
00701     invslope*=fWireTimetoCmCm*fWireTimetoCmCm;
00702     
00703     wireout=(ort_intercept - intercept)/(slope-invslope); 
00704     timeout=slope*wireout+intercept; 
00705   
00706     
00707     wireout/=fWiretoCm;
00708     timeout/=fTimetoCm;
00709     
00710     return 0;  
00711   }

Int_t larutil::GeometryUtilities::GetProjectedPoint ( const PxPoint *  p0,
const PxPoint *  p1,
PxPoint &  pN 
) const

Definition at line 823 of file GeometryUtilities.cc.

References larutil::Geometry::ChannelsIntersect(), detp, fNPlanes, fTimetoCm, geom, Get2DPointProjection(), larutil::Geometry::PlaneOriginVtx(), larutil::Geometry::PlaneWireToChannel(), and larutil::DetectorProperties::TriggerOffset().

00826   {
00827 
00828     //determine third plane number
00829     for(UInt_t i = 0; i < fNPlanes; ++i){
00830       if(i == p0->plane || i == p1->plane)
00831     continue;   
00832       pN.plane = i;
00833     }
00834   
00835     // Assuming there is no problem ( and we found the best pair that comes close in time )
00836     // we try to get the Y and Z coordinates for the start of the shower. 
00837     UInt_t chan1 = geom->PlaneWireToChannel(p0->plane,p0->w);
00838     UInt_t chan2 = geom->PlaneWireToChannel(p1->plane,p1->w);
00839  
00840     Double_t pos[3]={0.};
00841     geom->PlaneOriginVtx(p0->plane, pos);
00842  
00843     Double_t x=(p0->t - detp->TriggerOffset())*fTimetoCm+pos[0];
00844  
00845     Double_t y,z;
00846     if(! geom->ChannelsIntersect(chan1,chan2,y,z) )
00847       return -1;
00848  
00849  
00850     pos[1]=y;
00851     pos[2]=z;
00852     pos[0]=x;
00853     
00854     pN=Get2DPointProjection(pos, pN.plane);
00855        
00856     return 0;  
00857   }

Double_t larutil::GeometryUtilities::GetTimeTicks ( Double_t  x,
Int_t  plane 
) const

Definition at line 960 of file GeometryUtilities.cc.

References detp, fDriftVelocity, fTimeTick, geom, larutil::Geometry::PlaneOriginVtx(), and larutil::DetectorProperties::TriggerOffset().

00960                                                                        {
00961   
00962    
00963     
00964     Double_t pos[3];
00965     geom->PlaneOriginVtx(plane,pos);
00966     Double_t drifttick=(x/fDriftVelocity)*(1./fTimeTick);
00967     
00968     return drifttick-(pos[0]/fDriftVelocity)*(1./fTimeTick)+detp->TriggerOffset();  
00969     
00970      
00971    }

Int_t larutil::GeometryUtilities::GetYZ ( const PxPoint *  p0,
const PxPoint *  p1,
Double_t *  yz 
) const

Definition at line 861 of file GeometryUtilities.cc.

References larutil::Geometry::ChannelsIntersect(), geom, and larutil::Geometry::PlaneWireToChannel().

00864   {
00865     Double_t y,z;
00866   
00867     UInt_t chan1 = geom->PlaneWireToChannel(p0->plane, p0->w);
00868     UInt_t chan2 = geom->PlaneWireToChannel(p1->plane, p1->w);
00869 
00870     if(! geom->ChannelsIntersect(chan1,chan2,y,z) )
00871       return -1;
00872   
00873     yz[0]=y;
00874     yz[1]=z;
00875   
00876     return 0;
00877   }

Double_t larutil::GeometryUtilities::PitchInView ( UInt_t  plane,
Double_t  phi,
Double_t  theta 
) const

Todo:
switch to using new Geometry::WireAngleToVertical(geo::View_t)
Todo:
and Geometry::WirePitch(geo::View_t) methods

Definition at line 978 of file GeometryUtilities.cc.

References geom, GetDirectionCosines(), larutil::Geometry::PlaneToView(), larutil::Geometry::WireAngleToVertical(), and larutil::Geometry::WirePitch().

00981   {
00982     Double_t dirs[3] = {0.};
00983     GetDirectionCosines(phi,theta,dirs); 
00984     
00987     Double_t wirePitch   = 0.;
00988     Double_t angleToVert = 0.;
00989    
00990     wirePitch = geom->WirePitch(0,1,plane);
00991     angleToVert = geom->WireAngleToVertical(geom->PlaneToView(plane)) - 0.5*TMath::Pi();
00992          
00993     //(sin(angleToVert),std::cos(angleToVert)) is the direction perpendicular to wire
00994     //fDir.front() is the direction of the track at the beginning of its trajectory
00995     Double_t cosgamma = TMath::Abs(TMath::Sin(angleToVert)*dirs[1] + 
00996                       TMath::Cos(angleToVert)*dirs[2]);
00997    
00998     //   //std::cout << " ---- cosgamma: " << angleToVert*180/TMath::Pi() << " d's: " << dirs[1]
00999     //  << " " << dirs[2] << " ph,th " << phi << " " << theta << std::endl; 
01000     if(cosgamma < 1.e-5) 
01001       throw LArUtilException("cosgamma is basically 0, that can't be right");
01002     
01003     return wirePitch/cosgamma;
01004   }

std::vector< size_t > larutil::GeometryUtilities::PolyOverlap ( std::vector< const larutil::PxHit * >  ordered_hits,
std::vector< size_t >  candidate_polygon 
)

Definition at line 1310 of file GeometryUtilities.cc.

References Clockwise().

Referenced by SelectPolygonHitList().

01311                                                                  {
01312 
01313     //loop over edges
01314     for ( unsigned int i=0; i<(candidate_polygon.size()-1); i++){
01315       double Ax = ordered_hits.at(candidate_polygon.at(i))->w;
01316       double Ay = ordered_hits.at(candidate_polygon.at(i))->t;
01317       double Bx = ordered_hits.at(candidate_polygon.at(i+1))->w;
01318       double By = ordered_hits.at(candidate_polygon.at(i+1))->t;
01319       //loop over edges that have not been checked yet...
01320       //only ones furhter down in polygon
01321       for ( unsigned int j=i+2; j<(candidate_polygon.size()-1); j++){
01322     //avoid consecutive segments:
01323     if ( candidate_polygon.at(i) == candidate_polygon.at(j+1) )
01324       continue;
01325     else{
01326       double Cx = ordered_hits.at(candidate_polygon.at(j))->w;
01327       double Cy = ordered_hits.at(candidate_polygon.at(j))->t;
01328       double Dx = ordered_hits.at(candidate_polygon.at(j+1))->w;
01329       double Dy = ordered_hits.at(candidate_polygon.at(j+1))->t;
01330       
01331       if ( (Clockwise(Ax,Ay,Cx,Cy,Dx,Dy) != Clockwise(Bx,By,Cx,Cy,Dx,Dy))
01332            and (Clockwise(Ax,Ay,Bx,By,Cx,Cy) != Clockwise(Ax,Ay,Bx,By,Dx,Dy)) ){
01333         size_t tmp = candidate_polygon.at(i+1);
01334         candidate_polygon.at(i+1) = candidate_polygon.at(j);
01335         candidate_polygon.at(j) = tmp;
01336         //check that last element is still first (to close circle...)
01337         candidate_polygon.at(candidate_polygon.size()-1) = candidate_polygon.at(0);
01338         //swapped polygon...now do recursion to make sure
01339         return PolyOverlap( ordered_hits, candidate_polygon);
01340       }//if crossing
01341     }
01342       }//second loop
01343     }//first loop
01344     //std::cout << std::endl;
01345     return candidate_polygon;
01346   }

void larlight::larlight_base::print ( MSG::Level  level,
std::string  msg 
) const [inline, inherited]

message print out method

Definition at line 56 of file larlight_base.hh.

00057     {if(_verbosity.at(level)) Message::send(level,msg);};

void larlight::larlight_base::print ( MSG::Level  level,
std::string  where,
std::string  msg 
) const [inline, inherited]

message print out method

Definition at line 52 of file larlight_base.hh.

References larlight::larlight_base::_verbosity, and larlight::Message::send().

Referenced by CalculatePitch(), CalculatePitchPolar(), larlight::storage_manager::create_data_ptr(), larutil::LArProperties::DriftVelocity(), larutil::LArUtilBase::LoadData(), main(), larlight::storage_manager::prepare_tree(), and larlight::storage_manager::read_event().

00053     {if(_verbosity.at(level)) Message::send(level,where,msg);};

void larutil::GeometryUtilities::Reconfigure (  ) 

Definition at line 24 of file GeometryUtilities.cc.

References detp, larutil::LArProperties::DriftVelocity(), larutil::LArProperties::Efield(), fDriftVelocity, fNPlanes, fTimeTick, fTimetoCm, fWirePitch, fWireTimetoCmCm, fWiretoCm, geom, larutil::LArProperties::GetME(), larutil::DetectorProperties::GetME(), larutil::Geometry::GetME(), larp, larutil::Geometry::Nplanes(), larutil::Geometry::PlaneToView(), larutil::DetectorProperties::SamplingRate(), larutil::LArProperties::Temperature(), vertangle, larutil::Geometry::WireAngleToVertical(), and larutil::Geometry::WirePitch().

Referenced by GeometryUtilities(), and larutil::LArUtilManager::ReconfigureUtilities().

00025   {
00026     geom = (larutil::Geometry*)(larutil::Geometry::GetME());
00027     detp = (larutil::DetectorProperties*)(larutil::DetectorProperties::GetME());
00028     larp = (larutil::LArProperties*)(larutil::LArProperties::GetME());
00029     
00030     fNPlanes = geom->Nplanes();
00031     vertangle.resize(fNPlanes);
00032     for(UInt_t ip=0;ip<fNPlanes;ip++)
00033       vertangle[ip]=geom->WireAngleToVertical(geom->PlaneToView(ip)) - TMath::Pi()/2; // wire angle
00034         
00035     fWirePitch = geom->WirePitch(0,1,0);
00036     fTimeTick=detp->SamplingRate()/1000.; 
00037     fDriftVelocity=larp->DriftVelocity(larp->Efield(),larp->Temperature());
00038     
00039     fWiretoCm=fWirePitch;
00040     fTimetoCm=fTimeTick*fDriftVelocity;
00041     fWireTimetoCmCm=fTimetoCm/fWirePitch;
00042   }

void larutil::GeometryUtilities::SelectLocalHitlist ( const std::vector< larutil::PxHit > &  hitlist,
std::vector< const larutil::PxHit * > &  hitlistlocal,
larutil::PxPoint &  startHit,
Double_t &  linearlimit,
Double_t &  ortlimit,
Double_t &  lineslopetest,
larutil::PxHit &  averageHit 
)

Definition at line 1097 of file GeometryUtilities.cc.

References Get2DDistance(), and GetPointOnLine().

01104   {
01105 
01106     hitlistlocal.clear();
01107     double locintercept=startHit.t - startHit.w * lineslopetest;
01108     double timesum = 0;
01109     double wiresum = 0;
01110     for(size_t i=0; i<hitlist.size(); ++i) {
01111 
01112       larutil::PxPoint hitonline;
01113       
01114       GetPointOnLine( lineslopetest, locintercept, (const larutil::PxHit*)(&hitlist.at(i)), hitonline );
01115 
01116       //calculate linear distance from start point and orthogonal distance from axis
01117       Double_t lindist=Get2DDistance((const larutil::PxPoint*)(&hitonline),(const larutil::PxPoint*)(&startHit));
01118       Double_t ortdist=Get2DDistance((const larutil::PxPoint*)(&hitlist.at(i)),(const larutil::PxPoint*)(&hitonline));
01119       
01120       
01121       if(lindist<linearlimit && ortdist<ortlimit){
01122         hitlistlocal.push_back((const larutil::PxHit*)(&(hitlist.at(i))));
01123         timesum += hitlist.at(i).t;
01124         wiresum += hitlist.at(i).w;
01125       }
01126       
01127       
01128     }
01129 
01130     averageHit.plane = startHit.plane;
01131     if(hitlistlocal.size())
01132     {
01133       averageHit.w = wiresum/(double)hitlistlocal.size();
01134       averageHit.t = timesum/((double) hitlistlocal.size());
01135     }
01136   }

void larutil::GeometryUtilities::SelectLocalHitlist ( const std::vector< larlight::hit * > &  hitlist,
std::vector< UInt_t > &  hitlistlocal_index,
Double_t  wire_start,
Double_t  time_start,
Double_t  linearlimit,
Double_t  ortlimit,
Double_t  lineslopetest 
)

Definition at line 1032 of file GeometryUtilities.cc.

References Get2DDistance(), GetPlaneAndTPC(), and GetPointOnLine().

01037   {
01038     hitlistlocal_index.clear();
01039     Double_t locintercept=time_start-wire_start*lineslopetest;
01040     
01041     for(size_t i=0; i<hitlist.size(); ++i) {
01042 
01043       Double_t time = hitlist.at(i)->PeakTime();
01044       UInt_t plane,wire;
01045       GetPlaneAndTPC(hitlist.at(i),plane,wire);
01046       
01047       Double_t wonline=wire,tonline=time;
01048       //gser.GetPointOnLine(lineslopetest,lineinterctest,wire,time,wonline,tonline);
01049       GetPointOnLine(lineslopetest,locintercept,wire,time,wonline,tonline);
01050       
01051       //calculate linear distance from start point and orthogonal distance from axis
01052       Double_t lindist=Get2DDistance(wonline,tonline,wire_start,time_start);
01053       Double_t ortdist=Get2DDistance(wire,time,wonline,tonline);
01054       
01056       
01057       if(lindist<linearlimit && ortdist<ortlimit){ 
01058         hitlistlocal_index.push_back((UInt_t)i);      
01059         //std::cout << " w,t: " << wire << " " << time << " calc time: " << wire*lineslopetest + locintercept  << " ws,ts " << wonline << " "<< tonline <<" "<< lindist << " " << ortdist  << " plane: " << plane << std::endl;     
01060       }
01061       
01062     }
01063   }

void larutil::GeometryUtilities::SelectLocalHitlist ( const std::vector< larlight::hit * > &  hitlist,
std::vector< larlight::hit * > &  hitlistlocal_index,
Double_t  wire_start,
Double_t  time_start,
Double_t  linearlimit,
Double_t  ortlimit,
Double_t  lineslopetest 
)

Definition at line 1065 of file GeometryUtilities.cc.

References Get2DDistance(), GetPlaneAndTPC(), and GetPointOnLine().

01070   {
01071     hitlistlocal.clear();
01072     Double_t locintercept=time_start-wire_start*lineslopetest;
01073     
01074     for(size_t i=0; i<hitlist.size(); ++i) {
01075 
01076       Double_t time = hitlist.at(i)->PeakTime();
01077       UInt_t plane,wire;
01078       GetPlaneAndTPC(hitlist.at(i),plane,wire);
01079       
01080       Double_t wonline=wire,tonline=time;
01081       //gser.GetPointOnLine(lineslopetest,lineinterctest,wire,time,wonline,tonline);
01082       GetPointOnLine(lineslopetest,locintercept,wire,time,wonline,tonline);
01083       
01084       //calculate linear distance from start point and orthogonal distance from axis
01085       Double_t lindist=Get2DDistance(wonline,tonline,wire_start,time_start);
01086       Double_t ortdist=Get2DDistance(wire,time,wonline,tonline);
01087       
01089       
01090       if(lindist<linearlimit && ortdist<ortlimit){
01091         hitlistlocal.push_back(hitlist.at(i));
01092           //std::cout << " w,t: " << wire << " " << time << " calc time: " << wire*lineslopetest + locintercept  << " ws,ts " << wonline << " "<< tonline <<" "<< lindist << " " << ortdist  << " plane: " << plane << std::endl;
01093         }
01094     }
01095   }

void larutil::GeometryUtilities::SelectPolygonHitList ( const std::vector< larutil::PxHit > &  hitlist,
std::vector< const larutil::PxHit * > &  hitlistlocal 
)

Definition at line 1139 of file GeometryUtilities.cc.

References detp, fTimetoCm, fWiretoCm, geom, larutil::DetectorProperties::NumberTimeSamples(), larutil::Geometry::Nwires(), and PolyOverlap().

01141   {
01142     if(!(hitlist.size())) {
01143       throw LArUtilException("Provided empty hit list!");
01144       return;
01145     }
01146 
01147     hitlistlocal.clear();
01148     unsigned char plane = (*hitlist.begin()).plane;
01149 
01150     // Define subset of hits to define polygon
01151     std::map<double,const larutil::PxHit*> hitmap;
01152     double qtotal=0;
01153     for(auto const &h : hitlist){
01154       
01155       hitmap.insert(std::pair<double,const larutil::PxHit*>(h.charge,&h));
01156       qtotal += h.charge;
01157 
01158     }
01159     double qintegral=0;
01160     std::vector<const larutil::PxHit*> ordered_hits;
01161     ordered_hits.reserve(hitlist.size());
01162     for(auto hiter = hitmap.rbegin();
01163     qintegral < qtotal*0.95;
01164     ++hiter) {
01165 
01166       qintegral += (*hiter).first;
01167       ordered_hits.push_back((*hiter).second);
01168 
01169     }
01170 
01171     // Define container to hold found polygon corner PxHit index & distance
01172     std::vector<size_t> hit_index(8,0);
01173     std::vector<double> hit_distance(8,1e9);
01174     
01175     // Loop over hits and find corner points in the plane view
01176     // Also fill corner edge points
01177     std::vector<larutil::PxPoint> edges(4,PxPoint(plane,0,0));
01178     double wire_max = geom->Nwires(plane) * fWiretoCm;
01179     double time_max = detp->NumberTimeSamples() * fTimetoCm;
01180 
01181     for(size_t index = 0; index<ordered_hits.size(); ++index) {
01182 
01183       if(ordered_hits.at(index)->t < 0 ||
01184      ordered_hits.at(index)->w < 0 ||
01185      ordered_hits.at(index)->t > time_max ||
01186      ordered_hits.at(index)->w > wire_max ) {
01187 
01188     throw LArUtilException(Form("Invalid wire/time (%g,%g) ... range is (0=>%g,0=>%g)",
01189                     ordered_hits.at(index)->w,
01190                     ordered_hits.at(index)->t,
01191                     wire_max,
01192                     time_max)
01193                    );
01194     return;
01195       }
01196 
01197       double dist = 0;
01198       
01199       // Comparison w/ (Wire,0)
01200       dist = ordered_hits.at(index)->t;
01201       if(dist < hit_distance.at(1)) {
01202     hit_distance.at(1) = dist;
01203     hit_index.at(1) = index;
01204     edges.at(0).t = ordered_hits.at(index)->t;
01205     edges.at(1).t = ordered_hits.at(index)->t;
01206       }
01207 
01208       // Comparison w/ (WireMax,Time)
01209       dist = wire_max - ordered_hits.at(index)->w;
01210       if(dist < hit_distance.at(3)) {
01211     hit_distance.at(3) = dist;
01212     hit_index.at(3) = index;
01213     edges.at(1).w = ordered_hits.at(index)->w;
01214     edges.at(2).w = ordered_hits.at(index)->w;
01215       }
01216 
01217       // Comparison w/ (Wire,TimeMax)
01218       dist = time_max - ordered_hits.at(index)->t;
01219       if(dist < hit_distance.at(5)) {
01220     hit_distance.at(5) = dist;
01221     hit_index.at(5) = index;
01222     edges.at(2).t = ordered_hits.at(index)->t;
01223     edges.at(3).t = ordered_hits.at(index)->t;
01224       }
01225 
01226       // Comparison w/ (0,Time)
01227       dist = ordered_hits.at(index)->w;
01228       if(dist < hit_distance.at(7)) {
01229     hit_distance.at(7) = dist;
01230     hit_index.at(7) = index;
01231     edges.at(0).w = ordered_hits.at(index)->w;
01232     edges.at(3).w = ordered_hits.at(index)->w;
01233       }
01234     }
01235     /*
01236     std::cout
01237       << Form("Corner points: (%g,%g) (%g,%g) (%g,%g) (%g,%g)",
01238           edges.at(0).w, edges.at(0).t,
01239           edges.at(1).w, edges.at(1).t,
01240           edges.at(2).w, edges.at(2).t,
01241           edges.at(3).w, edges.at(3).t)
01242       <<std::endl;
01243     */
01244     for(size_t index = 0; index<ordered_hits.size(); ++index) {
01245 
01246       double dist = 0;
01247       // Comparison w/ (0,0)
01248       dist = pow((ordered_hits.at(index)->t - edges.at(0).t),2) + pow((ordered_hits.at(index)->w - edges.at(0).w),2);
01249       if(dist < hit_distance.at(0)) {
01250     hit_distance.at(0) = dist;
01251     hit_index.at(0) = index;
01252       }
01253       
01254       // Comparison w/ (WireMax,0)
01255       dist = pow((ordered_hits.at(index)->t - edges.at(1).t),2) + pow((ordered_hits.at(index)->w - edges.at(1).w),2);
01256       if(dist < hit_distance.at(2)) {
01257     hit_distance.at(2) = dist;
01258     hit_index.at(2) = index;
01259       }
01260 
01261       // Comparison w/ (WireMax,TimeMax)
01262       dist = pow((ordered_hits.at(index)->t - edges.at(2).t),2) + pow((ordered_hits.at(index)->w - edges.at(2).w),2);
01263       if(dist < hit_distance.at(4)) {
01264     hit_distance.at(4) = dist;
01265     hit_index.at(4) = index;
01266       }
01267 
01268       // Comparison w/ (0,TimeMax)
01269       dist = pow((ordered_hits.at(index)->t - edges.at(3).t),2) + pow((ordered_hits.at(index)->w - edges.at(3).w),2);
01270       if(dist < hit_distance.at(6)) {
01271     hit_distance.at(6) = dist;
01272     hit_index.at(6) = index;
01273       }
01274 
01275     }
01276 
01277     // Loop over the resulting hit indexes and append unique hits to define the polygon to the return hit list
01278     std::set<size_t> unique_index;
01279     std::vector<size_t> candidate_polygon;
01280     candidate_polygon.reserve(9);
01281     //    std::cout << "Original polygon: " << std::endl;
01282     for(auto &index : hit_index) {
01283       
01284       if(unique_index.find(index) == unique_index.end()) {
01285     //  hitlistlocal.push_back((const larutil::PxHit*)(ordered_hits.at(index)));
01286     //std::cout << "(" << ordered_hits.at(index)->w << ", " << ordered_hits.at(index)->t << ")" << std::endl;
01287     unique_index.insert(index);
01288     candidate_polygon.push_back(index);
01289       }
01290     }
01291     for (auto &index: hit_index){
01292       candidate_polygon.push_back(index);
01293       break;
01294     }
01295 
01296     if(unique_index.size()>8) throw LArUtilException("Size of the polygon > 8!");    
01297 
01298     //Untangle Polygon
01299     candidate_polygon = PolyOverlap( ordered_hits, candidate_polygon);
01300     
01301     hitlistlocal.clear();
01302     for( unsigned int i=0; i<(candidate_polygon.size()-1); i++){
01303       hitlistlocal.push_back((const larutil::PxHit*)(ordered_hits.at(candidate_polygon.at(i))));
01304     }
01305     //check that polygon does not have more than 8 sides
01306     if(unique_index.size()>8) throw LArUtilException("Size of the polygon > 8!");    
01307   }

void larlight::larlight_base::set_verbosity ( MSG::Level  level  )  [virtual, inherited]

Setter for the verbosity level.

Reimplemented in larlight::ana_processor.

Definition at line 8 of file larlight_base.cc.

References larlight::larlight_base::_verbosity, larlight::larlight_base::_verbosity_level, larlight::MSG::DEBUG, larlight::MSG::ERROR, larlight::MSG::INFO, larlight::MSG::MSG_TYPE_MAX, larlight::MSG::NORMAL, and larlight::MSG::WARNING.

Referenced by larlight::larlight_base::larlight_base(), main(), and larlight::ana_processor::set_verbosity().

00010   {
00011     
00012     _verbosity_level=level;
00013     
00014     for(size_t i=(size_t)(MSG::DEBUG); i<(size_t)(MSG::MSG_TYPE_MAX); ++i)
00015       _verbosity[i]=false;
00016     
00017     switch(level){
00018     case MSG::DEBUG:
00019       _verbosity[MSG::DEBUG]=true;
00020     case MSG::INFO:
00021       _verbosity[MSG::INFO]=true;
00022     case MSG::NORMAL:
00023       _verbosity[MSG::NORMAL]=true;
00024     case MSG::WARNING:
00025       _verbosity[MSG::WARNING]=true;
00026     case MSG::ERROR:
00027       _verbosity[MSG::ERROR]=true;
00028     case MSG::MSG_TYPE_MAX:
00029       break;
00030     }
00031     
00032   }

Double_t larutil::GeometryUtilities::TimeToCm (  )  const [inline]

Definition at line 228 of file GeometryUtilities.hh.

References fTimetoCm.

00228 {return fTimetoCm;}

Double_t larutil::GeometryUtilities::WireTimeToCmCm (  )  const [inline]

Definition at line 230 of file GeometryUtilities.hh.

References fWireTimetoCmCm.

00230 {return fWireTimetoCmCm;}

Double_t larutil::GeometryUtilities::WireToCm (  )  const [inline]

Definition at line 229 of file GeometryUtilities.hh.

References fWiretoCm.

00229 {return fWiretoCm;}

Member Data Documentation

char larlight::larlight_base::_buf[200] [protected, inherited]

char buffer for message manipulation

Definition at line 57 of file larlight_base.hh.

Referenced by larlight::storage_manager::open(), larlight::storage_manager::prepare_tree(), and larlight::ana_processor::run().

GeometryUtilities * larutil::GeometryUtilities::_me = 0 [static, private]

Definition at line 36 of file GeometryUtilities.hh.

Referenced by GetME().

std::string larlight::larlight_base::_name [protected, inherited]

class name holder

Definition at line 64 of file larlight_base.hh.

Referenced by larlight::ana_base::ana_base(), larlight::ana_processor::ana_processor(), larlight::larlight_base::class_name(), larutil::DetectorProperties::DetectorProperties(), larutil::Geometry::Geometry(), GeometryUtilities(), larutil::LArProperties::LArProperties(), larutil::LArUtilBase::LArUtilBase(), larutil::LArUtilBase::LoadData(), and larlight::storage_manager::storage_manager().

std::vector<bool> larlight::larlight_base::_verbosity [protected, inherited]

holder for enabled message levels

Definition at line 62 of file larlight_base.hh.

Referenced by larlight::storage_manager::close(), larlight::ana_processor::finalize(), larlight::ana_processor::initialize(), larlight::storage_manager::open(), larlight::storage_manager::prepare_tree(), larlight::larlight_base::print(), larlight::storage_manager::reset(), larlight::ana_processor::reset(), larlight::ana_processor::run(), and larlight::larlight_base::set_verbosity().

MSG::Level larlight::larlight_base::_verbosity_level [protected, inherited]

holder for specified verbosity level

Definition at line 63 of file larlight_base.hh.

Referenced by larlight::larlight_base::get_verbosity(), larlight::ana_processor::initialize(), and larlight::larlight_base::set_verbosity().

larutil::DetectorProperties* larutil::GeometryUtilities::detp [private]

Definition at line 235 of file GeometryUtilities.hh.

Referenced by Get2DPointProjection(), GetProjectedPoint(), GetTimeTicks(), Reconfigure(), and SelectPolygonHitList().

Double_t larutil::GeometryUtilities::fDriftVelocity [private]

Definition at line 241 of file GeometryUtilities.hh.

Referenced by Get2DPointProjection(), GetTimeTicks(), and Reconfigure().

UInt_t larutil::GeometryUtilities::fNPlanes [private]

Definition at line 242 of file GeometryUtilities.hh.

Referenced by GetProjectedPoint(), and Reconfigure().

Double_t larutil::GeometryUtilities::fTimeTick [private]

Definition at line 240 of file GeometryUtilities.hh.

Referenced by Get2DPointProjection(), GetTimeTicks(), and Reconfigure().

Double_t larutil::GeometryUtilities::fTimetoCm [private]

Definition at line 244 of file GeometryUtilities.hh.

Referenced by Get2Dangle(), Get2DDistance(), Get2Dslope(), GetPointOnLineWSlopes(), GetProjectedPoint(), Reconfigure(), SelectPolygonHitList(), and TimeToCm().

Double_t larutil::GeometryUtilities::fWirePitch [private]

Definition at line 239 of file GeometryUtilities.hh.

Referenced by Reconfigure().

Double_t larutil::GeometryUtilities::fWireTimetoCmCm [private]

Definition at line 245 of file GeometryUtilities.hh.

Referenced by Get2Dslope(), GetPointOnLineWSlopes(), Reconfigure(), and WireTimeToCmCm().

Double_t larutil::GeometryUtilities::fWiretoCm [private]

Definition at line 243 of file GeometryUtilities.hh.

Referenced by Get2Dangle(), Get2DDistance(), Get2DPitchDistance(), Get2DPitchDistanceWSlope(), Get2DPointProjectionCM(), Get2Dslope(), GetPointOnLineWSlopes(), Reconfigure(), SelectPolygonHitList(), and WireToCm().

larutil::Geometry* larutil::GeometryUtilities::geom [private]

Definition at line 234 of file GeometryUtilities.hh.

Referenced by CalculatePitch(), CalculatePitchPolar(), Get2DangleFrom3D(), Get2DPointProjection(), Get2DPointProjectionCM(), GetPlaneAndTPC(), GetProjectedPoint(), GetTimeTicks(), GetYZ(), PitchInView(), Reconfigure(), and SelectPolygonHitList().

larutil::LArProperties* larutil::GeometryUtilities::larp [private]

Definition at line 236 of file GeometryUtilities.hh.

Referenced by Reconfigure().

std::vector< Double_t > larutil::GeometryUtilities::vertangle [private]

Definition at line 238 of file GeometryUtilities.hh.

Referenced by Get3DaxisN(), Get3DSpecialCaseTheta(), and Reconfigure().

The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Enumerations Enumerator
Generated on 3 Jun 2014 for MyProject by  doxygen 1.6.1