CMS 3D CMS Logo

CachedTrajectory Class Reference

#include <TrackingTools/TrackAssociator/interface/CachedTrajectory.h>

List of all members.

Public Types

enum  TrajectorType { IpToEcal, IpToHcal, IpToHO, FullTrajectory }
enum  WideTrajectoryType { Ecal, Hcal, HO }

Public Member Functions

 CachedTrajectory ()
void findEcalTrajectory (const FiducialVolume &)
void findHcalTrajectory (const FiducialVolume &)
void findHOTrajectory (const FiducialVolume &)
const std::vector
< SteppingHelixStateInfo > & 
getEcalTrajectory ()
const std::vector
< SteppingHelixStateInfo > & 
getHcalTrajectory ()
const std::vector
< SteppingHelixStateInfo > & 
getHOTrajectory ()
SteppingHelixStateInfo getInnerState ()
SteppingHelixStateInfo getOuterState ()
float getPropagationStep () const
SteppingHelixStateInfo getStateAtEcal ()
SteppingHelixStateInfo getStateAtHcal ()
SteppingHelixStateInfo getStateAtHO ()
void getTrajectory (std::vector< SteppingHelixStateInfo > &, const FiducialVolume &, int steps=4)
 get a set of points representing the trajectory between two cylinders of radius R1 and R2 and length L1 and L2.
std::vector< GlobalPoint > * getWideTrajectory (const std::vector< SteppingHelixStateInfo > &, WideTrajectoryType)
TrajectoryStateOnSurface propagate (const Plane *plane)
 get fast to a given DetId surface using cached trajectory
bool propagateAll (const SteppingHelixStateInfo &initialState)
 propagate through the whole detector, returns true if successful
void propagateForward (SteppingHelixStateInfo &state, float distance)
void reset_trajectory ()
void setMaxDetectorLength (float l=2200.)
void setMaxDetectorRadius (float r=800.)
void setMaxHOLength (float l=2200.)
void setMaxHORadius (float r=800.)
void setMinDetectorLength (float l=0.)
void setMinDetectorRadius (float r=0.)
void setPropagationStep (float s=20.)
void setPropagator (const Propagator *ptr)
void setStateAtIP (const SteppingHelixStateInfo &state)
std::pair< float, float > trajectoryDelta (TrajectorType)
 calculate trajectory change (Theta,Phi) delta = final - original

Protected Member Functions

std::pair< float, float > delta (const double &theta1, const double &theta2, const double &phi1, const double &phi2)
float distance (const Plane *plane, int index)

Static Protected Member Functions

static int sign (float number)

Protected Attributes

std::vector
< SteppingHelixStateInfo
ecalTrajectory_
std::deque
< SteppingHelixStateInfo
fullTrajectory_
bool fullTrajectoryFilled_
std::vector
< SteppingHelixStateInfo
hcalTrajectory_
float HOmaxRho_
float HOmaxZ_
std::vector
< SteppingHelixStateInfo
hoTrajectory_
float maxRho_
float maxZ_
float minRho_
float minZ_
const Propagatorpropagator_
SteppingHelixStateInfo stateAtIP_
float step_
std::vector< GlobalPointwideEcalTrajectory_
std::vector< GlobalPointwideHcalTrajectory_
std::vector< GlobalPointwideHOTrajectory_


Detailed Description

Definition at line 34 of file CachedTrajectory.h.


Member Enumeration Documentation

enum CachedTrajectory::TrajectorType

Enumerator:
IpToEcal 
IpToHcal 
IpToHO 
FullTrajectory 

Definition at line 37 of file CachedTrajectory.h.

enum CachedTrajectory::WideTrajectoryType

Enumerator:
Ecal 
Hcal 
HO 

Definition at line 38 of file CachedTrajectory.h.

00038 { Ecal, Hcal, HO };


Constructor & Destructor Documentation

CachedTrajectory::CachedTrajectory (  ) 

Definition at line 19 of file CachedTrajectory.cc.

References reset_trajectory(), setMaxDetectorLength(), setMaxDetectorRadius(), setMinDetectorLength(), setMinDetectorRadius(), and setPropagationStep().


Member Function Documentation

std::pair< float, float > CachedTrajectory::delta ( const double &  theta1,
const double &  theta2,
const double &  phi1,
const double &  phi2 
) [protected]

Definition at line 279 of file CachedTrajectory.cc.

References HLT_VtxMuL3::result.

Referenced by trajectoryDelta().

00283 {
00284    std::pair<float,float> result(theta2 - theta1, phi2 - phi1 );
00285    // this won't work for loopers, since deltaPhi cannot be larger than Pi.
00286    if ( fabs(result.second) > 2*M_PI-fabs(result.second) ) {
00287       if (result.second>0) 
00288         result.second -= 2*M_PI;
00289       else
00290         result.second += 2*M_PI;
00291    }
00292    return result;
00293 }

float CachedTrajectory::distance ( const Plane plane,
int  index 
) [inline, protected]

Definition at line 113 of file CachedTrajectory.h.

References fullTrajectory_, and Plane::localZ().

Referenced by propagate().

00113                                                  {
00114       if (index<0 || fullTrajectory_.empty() || uint(index) >= fullTrajectory_.size()) return 0;
00115       return plane->localZ(fullTrajectory_[index].position());
00116    }

void CachedTrajectory::findEcalTrajectory ( const FiducialVolume volume  ) 

Definition at line 376 of file CachedTrajectory.cc.

References ecalTrajectory_, getTrajectory(), and LogTrace.

Referenced by TrackDetectorAssociator::associate().

00376                                                                         {
00377    LogTrace("TrackAssociator") << "getting trajectory in ECAL";
00378    getTrajectory(ecalTrajectory_, volume, 4 );
00379    LogTrace("TrackAssociator") << "# of points in ECAL trajectory:" << ecalTrajectory_.size();
00380 }

void CachedTrajectory::findHcalTrajectory ( const FiducialVolume volume  ) 

Definition at line 386 of file CachedTrajectory.cc.

References getTrajectory(), hcalTrajectory_, and LogTrace.

Referenced by TrackDetectorAssociator::associate().

00386                                                                         {
00387    LogTrace("TrackAssociator") << "getting trajectory in HCAL";
00388    getTrajectory(hcalTrajectory_, volume, 20 ); // more steps to account for different depth
00389    LogTrace("TrackAssociator") << "# of points in HCAL trajectory:" << hcalTrajectory_.size();
00390 }

void CachedTrajectory::findHOTrajectory ( const FiducialVolume volume  ) 

Definition at line 396 of file CachedTrajectory.cc.

References getTrajectory(), hoTrajectory_, and LogTrace.

Referenced by TrackDetectorAssociator::associate().

00396                                                                       {
00397    LogTrace("TrackAssociator") << "getting trajectory in HO";
00398    getTrajectory(hoTrajectory_, volume, 2 );
00399    LogTrace("TrackAssociator") << "# of points in HO trajectory:" << hoTrajectory_.size();
00400 }

const std::vector< SteppingHelixStateInfo > & CachedTrajectory::getEcalTrajectory (  ) 

Definition at line 382 of file CachedTrajectory.cc.

References ecalTrajectory_.

Referenced by TrackDetectorAssociator::fillCaloTowers(), and TrackDetectorAssociator::fillEcal().

00382                                                                              {
00383    return ecalTrajectory_;
00384 }

const std::vector< SteppingHelixStateInfo > & CachedTrajectory::getHcalTrajectory (  ) 

Definition at line 392 of file CachedTrajectory.cc.

References hcalTrajectory_.

Referenced by TrackDetectorAssociator::fillCaloTowers(), and TrackDetectorAssociator::fillHcal().

00392                                                                              {
00393    return hcalTrajectory_;
00394 }

const std::vector< SteppingHelixStateInfo > & CachedTrajectory::getHOTrajectory (  ) 

Definition at line 402 of file CachedTrajectory.cc.

References hoTrajectory_.

Referenced by TrackDetectorAssociator::fillHO().

00402                                                                            {
00403    return hoTrajectory_;
00404 }

SteppingHelixStateInfo CachedTrajectory::getInnerState (  ) 

Definition at line 526 of file CachedTrajectory.cc.

References fullTrajectory_.

00526                                                        {
00527   if(fullTrajectory_.empty() )
00528     return SteppingHelixStateInfo();
00529   else
00530     return fullTrajectory_.front();
00531 }

SteppingHelixStateInfo CachedTrajectory::getOuterState (  ) 

Definition at line 534 of file CachedTrajectory.cc.

References fullTrajectory_.

Referenced by TrackDetectorAssociator::getTAMuonChamberMatches().

00534                                                        {
00535   if(fullTrajectory_.empty() )
00536     return SteppingHelixStateInfo();
00537   else
00538     return fullTrajectory_.back();
00539 }

float CachedTrajectory::getPropagationStep (  )  const [inline]

Definition at line 96 of file CachedTrajectory.h.

References step_.

Referenced by TrackDetectorAssociator::associate().

00096 { return step_;}

SteppingHelixStateInfo CachedTrajectory::getStateAtEcal (  ) 

Definition at line 501 of file CachedTrajectory.cc.

References ecalTrajectory_.

Referenced by TrackDetectorAssociator::associate().

00502 {
00503    if ( ecalTrajectory_.empty() )
00504      return SteppingHelixStateInfo();
00505    else 
00506      return ecalTrajectory_.front();
00507 }

SteppingHelixStateInfo CachedTrajectory::getStateAtHcal (  ) 

Definition at line 509 of file CachedTrajectory.cc.

References hcalTrajectory_.

Referenced by TrackDetectorAssociator::associate(), and TrackDetectorAssociator::getTAMuonChamberMatches().

00510 {
00511    if ( hcalTrajectory_.empty() )
00512      return SteppingHelixStateInfo();
00513    else 
00514      return hcalTrajectory_.front();
00515 }

SteppingHelixStateInfo CachedTrajectory::getStateAtHO (  ) 

Definition at line 517 of file CachedTrajectory.cc.

References hoTrajectory_.

Referenced by TrackDetectorAssociator::associate().

00518 {
00519    if ( hoTrajectory_.empty() )
00520      return SteppingHelixStateInfo();
00521    else 
00522      return hoTrajectory_.front();
00523 }

void CachedTrajectory::getTrajectory ( std::vector< SteppingHelixStateInfo > &  trajectory,
const FiducialVolume volume,
int  steps = 4 
)

get a set of points representing the trajectory between two cylinders of radius R1 and R2 and length L1 and L2.

Parameter steps defines maximal number of steps in the detector.

Definition at line 295 of file CachedTrajectory.cc.

References Exception, fullTrajectory_, fullTrajectoryFilled_, i, FiducialVolume::isValid(), SteppingHelixStateInfo::isValid(), LogTrace, max, FiducialVolume::maxR(), FiducialVolume::maxZ(), FiducialVolume::minR(), FiducialVolume::minZ(), SteppingHelixStateInfo::momentum(), PV3DBase< T, PVType, FrameType >::perp(), SteppingHelixStateInfo::position(), propagateForward(), cmsRelvalreportInput::step, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), PV3DBase< T, PVType, FrameType >::z(), and z.

Referenced by findEcalTrajectory(), findHcalTrajectory(), and findHOTrajectory().

00298 {
00299    if ( ! fullTrajectoryFilled_ ) throw cms::Exception("FatalError") << "trajectory is not defined yet. Please use propagateAll first.";
00300    if ( fullTrajectory_.empty() ) {
00301       LogTrace("TrackAssociator") << "Trajectory is empty. Move on";
00302       return;
00303    }
00304    if ( ! volume.isValid() ) {
00305       LogTrace("TrackAssociator") << "no trajectory is expected to be found since the fiducial volume is not valid";
00306       return;
00307    }
00308    double step = std::max(volume.maxR()-volume.minR(),volume.maxZ()-volume.minZ())/steps;
00309    
00310    int closestPointOnLeft = -1;
00311    
00312    // check whether the trajectory crossed the region
00313    if ( ! 
00314         ( ( fullTrajectory_.front().position().perp() < volume.maxR() && fabs(fullTrajectory_.front().position().z()) < volume.maxZ() ) &&
00315           ( fullTrajectory_.back().position().perp() > volume.minR() || fabs(fullTrajectory_.back().position().z()) > volume.minZ() ) ))
00316      {
00317         LogTrace("TrackAssociator") << "Track didn't cross the region (R1,R2,L1,L2): " << volume.minR() << ", " << volume.maxR() <<
00318           ", " << volume.minZ() << ", " << volume.maxZ();
00319         return;
00320      }
00321    
00322    // get distance along momentum to the surface.
00323    
00324    // the following code can be made faster, but it'll hardly be a significant improvement
00325    // simplifications:
00326    //   1) direct loop over stored trajectory points instead of some sort 
00327    //      of fast root search (Newton method)
00328    //   2) propagate from the closest point outside the region with the 
00329    //      requested step ignoring stored trajectory points.
00330    for(uint i=0; i<fullTrajectory_.size(); i++) {
00331       // LogTrace("TrackAssociator") << "Trajectory info (i,perp,r1,r2,z,z1,z2): " << i << ", " << fullTrajectory_[i].position().perp() <<
00332       //        ", " << volume.minR() << ", " << volume.maxR() << ", " << fullTrajectory_[i].position().z() << ", " << volume.minZ() << ", " << 
00333       //        volume.maxZ() << ", " << closestPointOnLeft;
00334       if ( fullTrajectory_[i].position().perp()-volume.minR() > 0  || fabs(fullTrajectory_[i].position().z()) - volume.minZ() >0 )
00335         {
00336            if (i>0) 
00337              closestPointOnLeft = i - 1;
00338            else
00339              closestPointOnLeft = 0;
00340            break;
00341         }
00342    }
00343    if (closestPointOnLeft == -1) throw cms::Exception("FatalError") << "This shouls never happen - internal logic error";
00344    
00345    SteppingHelixStateInfo currentState(fullTrajectory_[closestPointOnLeft]);
00346    if ( currentState.position().x()*currentState.momentum().x() +
00347         currentState.position().y()*currentState.momentum().y() +
00348         currentState.position().z()*currentState.momentum().z() < 0 )
00349      step = -step;
00350    
00351    while (currentState.position().perp() < volume.maxR() && fabs(currentState.position().z()) < volume.maxZ() )
00352      {
00353         propagateForward(currentState,step);
00354         if (! currentState.isValid() ) {
00355            LogTrace("TrackAssociator") << "Failed to propagate the track; moving on\n";
00356            break;
00357         }
00358         // LogTrace("TrackAssociator") << "New state (perp, z): " << currentState.position().perp() << ", " << currentState.position().z();
00359         if ( ( currentState.position().perp() < volume.maxR() && fabs(currentState.position().z()) < volume.maxZ() ) &&
00360              ( currentState.position().perp()-volume.minR() > 0  || fabs(currentState.position().z()) - volume.minZ() >0 ) )
00361           trajectory.push_back(currentState);
00362      }
00363 }

std::vector< GlobalPoint > * CachedTrajectory::getWideTrajectory ( const std::vector< SteppingHelixStateInfo > &  states,
WideTrajectoryType  wideTrajectoryType 
)

Definition at line 407 of file CachedTrajectory.cc.

References funct::cos(), Ecal, Hcal, HO, index, LogTrace, M_PI_2, phi, funct::pow(), dttmaxenums::r32, LocalError::rotate(), LocalError::scale(), funct::sin(), funct::sqrt(), state, target, Surface::toGlobal(), wideEcalTrajectory_, wideHcalTrajectory_, wideHOTrajectory_, LocalError::xx(), LocalError::xy(), and LocalError::yy().

Referenced by TrackDetectorAssociator::fillEcal(), TrackDetectorAssociator::fillHcal(), and TrackDetectorAssociator::fillHO().

00408                                                                            {
00409    std::vector<GlobalPoint>* wideTrajectory = 0;
00410    switch (wideTrajectoryType) {
00411     case Ecal:
00412        LogTrace("TrackAssociator") << "Filling ellipses in Ecal trajectory";
00413        wideTrajectory = &wideEcalTrajectory_;
00414        break;
00415     case Hcal:
00416        LogTrace("TrackAssociator") << "Filling ellipses in Hcal trajectory";
00417       wideTrajectory = &wideHcalTrajectory_;
00418       break;
00419     case HO:
00420        LogTrace("TrackAssociator") << "Filling ellipses in HO trajectory";
00421        wideTrajectory = &wideHOTrajectory_;
00422        break;
00423    }
00424    if(!wideTrajectory) return 0;
00425 
00426    for(std::vector<SteppingHelixStateInfo>::const_iterator state= states.begin();
00427        state != states.end(); state++) {
00428       // defined a normal plane wrt the particle trajectory direction
00429       // let's hope that I computed the rotation matrix correctly.
00430       GlobalVector vector(state->momentum().unit());
00431       float r21 = 0;
00432       float r22 = vector.z()/sqrt(1-pow(vector.x(),2));
00433       float r23 = -vector.y()/sqrt(1-pow(vector.x(),2));
00434       float r31 = vector.x();
00435       float r32 = vector.y();
00436       float r33 = vector.z();
00437       float r11 = r22*r33-r23*r32;
00438       float r12 = r23*r31;
00439       float r13 = -r22*r31;
00440    
00441       Plane::RotationType rotation(r11, r12, r13,
00442                                    r21, r22, r23,
00443                                    r31, r32, r33);
00444       Plane* target = new Plane(state->position(), rotation);
00445 
00446       TrajectoryStateOnSurface tsos = state->getStateOnSurface(*target);
00447 
00448       if (!tsos.isValid()) {
00449          LogTrace("TrackAssociator") << "[getWideTrajectory] TSOS not valid";
00450          continue;
00451       }
00452       if (!tsos.hasError()) {
00453          LogTrace("TrackAssociator") << "[getWideTrajectory] TSOS does not have Errors";
00454          continue;
00455       }
00456       LocalError localErr = tsos.localError().positionError();
00457       localErr.scale(2); // get the 2 sigma ellipse
00458       float xx = localErr.xx();
00459       float xy = localErr.xy();
00460       float yy = localErr.yy();
00461 
00462       float denom = yy - xx;
00463       float phi = 0., phi_temp=0.;
00464       if(xy == 0 && denom==0) phi = M_PI_4;
00465       else phi = 0.5 * atan2(2.*xy,denom); // angle of MAJOR axis
00466       phi_temp = phi;
00467       // Unrotate the error ellipse to get the semimajor and minor axes. Then place points on
00468       // the endpoints of semiminor an seminajor axes on original(rotated) error ellipse.
00469       LocalError rotErr = localErr.rotate(-phi); // xy covariance of rotErr should be zero
00470       float semi1 = sqrt(rotErr.xx());
00471       float semi2 = sqrt(rotErr.yy());
00472       
00473       // Just use one point if the ellipse is small
00474       // if(semi1 < 0.1 && semi2 < 0.1) {
00475       //   LogTrace("TrackAssociator") << "[getWideTrajectory] Error ellipse is small, using one trajectory point";
00476       //   wideTrajectory->push_back(state->position());
00477       //   continue;
00478       // }
00479 
00480       Local2DPoint bounds[4];
00481       bounds[0] = Local2DPoint(semi1*cos(phi),         semi1*sin(phi));
00482       bounds[1] = Local2DPoint(semi1*cos(phi+M_PI),    semi1*sin(phi+M_PI));
00483       phi += M_PI_2; // add pi/2 for the semi2 axis
00484       bounds[2] = Local2DPoint(semi2*cos(phi),         semi2*sin(phi));
00485       bounds[3] = Local2DPoint(semi2*cos(phi+M_PI),    semi2*sin(phi+M_PI));
00486 
00487       // LogTrace("TrackAssociator") << "Axes " << semi1 <<","<< semi2 <<"   phi "<< phi;
00488       // LogTrace("TrackAssociator") << "Local error ellipse: " << bounds[0] << bounds[1] << bounds[2] << bounds[3];
00489 
00490       wideTrajectory->push_back(state->position());
00491       for(int index=0; index<4; ++index)
00492          wideTrajectory->push_back(target->toGlobal(bounds[index]));
00493 
00494       //LogTrace("TrackAssociator") <<"Global error ellipse: (" << target->toGlobal(bounds[0]) <<","<< target->toGlobal(bounds[1])
00495       //         <<","<< target->toGlobal(bounds[2]) <<","<< target->toGlobal(bounds[3]) <<","<<state->position() <<")";
00496    }
00497 
00498    return wideTrajectory;
00499 }

TrajectoryStateOnSurface CachedTrajectory::propagate ( const Plane plane  ) 

get fast to a given DetId surface using cached trajectory

Definition at line 142 of file CachedTrajectory.cc.

References cms::Exception::category(), dist(), distance(), cms::Exception::explainSelf(), fullTrajectory_, SteppingHelixStateInfo::getStateOnSurface(), int, LogTrace, SteppingHelixStateInfo::momentum(), GloballyPositioned< T >::position(), Propagator::propagate(), propagator_, sign(), state, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by TrackDetectorAssociator::getTAMuonChamberMatches().

00143 {
00144    // TimerStack timers(TimerStack::Disableable);
00145    // timers.benchmark("CachedTrajectory::propagate::benchmark");
00146    // timers.push("CachedTrajectory::propagate",TimerStack::FastMonitoring);
00147    // timers.push("CachedTrajectory::propagate::findClosestPoint",TimerStack::FastMonitoring);
00148 
00149    // Assume that all points along the trajectory are equally spread out.
00150    // For simplication assume that the trajectory is just a straight
00151    // line and find a point closest to the target plane. Propagate to
00152    // the plane from the point.
00153    
00154    const float matchingDistance = 1;
00155    // find the closest point to the plane
00156    int leftIndex = 0;
00157    int rightIndex = fullTrajectory_.size()-1;
00158    int closestPointOnLeft = 0;
00159    
00160    // check whether the trajectory crossed the plane (signs should be different)
00161    if ( sign( distance(plane, leftIndex) ) * sign( distance(plane, rightIndex) ) != -1 ) {
00162       LogTrace("TrackAssociator") << "Track didn't cross the plane:\n\tleft distance: "<<distance(plane, leftIndex)
00163         <<"\n\tright distance: " << distance(plane, rightIndex);
00164      return TrajectoryStateOnSurface();
00165    }
00166    
00167    while (leftIndex + 1 < rightIndex) {
00168       closestPointOnLeft = int((leftIndex+rightIndex)/2);
00169       float dist = distance(plane,closestPointOnLeft);
00170       // LogTrace("TrackAssociator") << "Closest point on left: " << closestPointOnLeft << "\n"
00171       //    << "Distance to the plane: " << dist;
00172       if (fabs(dist)<matchingDistance) {
00173          // found close match, verify that we are on the left side
00174          if (closestPointOnLeft>0 && sign( distance(plane, closestPointOnLeft-1) ) * dist == -1)
00175            closestPointOnLeft--;
00176          break; 
00177       }
00178       
00179       // check where is the plane
00180       if (sign( distance(plane, leftIndex) * dist ) == -1)
00181         rightIndex = closestPointOnLeft;
00182       else
00183         leftIndex = closestPointOnLeft;
00184       
00185       // LogTrace("TrackAssociator") << "Distance on left: " << distance(plane, leftIndex) << "\n"
00186       //        << "Distance to closest point: " <<  distance(plane, closestPointOnLeft) << "\n"
00187       //        << "Left index: " << leftIndex << "\n"
00188       //        << "Right index: " << rightIndex;
00189       
00190    }
00191       LogTrace("TrackAssociator") << "closestPointOnLeft: " << closestPointOnLeft 
00192         << "\n\ttrajectory point (z,R,eta,phi): " 
00193         << fullTrajectory_[closestPointOnLeft].position().z() << ", "
00194         << fullTrajectory_[closestPointOnLeft].position().perp() << " , "       
00195         << fullTrajectory_[closestPointOnLeft].position().eta() << " , " 
00196         << fullTrajectory_[closestPointOnLeft].position().phi()
00197         << "\n\tplane center (z,R,eta,phi): " 
00198         << plane->position().z() << ", "
00199         << plane->position().perp() << " , "    
00200         << plane->position().eta() << " , " 
00201         << plane->position().phi();
00202      
00203    // propagate to the plane
00204    // timers.pop_and_push("CachedTrajectory::propagate::localPropagation",TimerStack::FastMonitoring);
00205    if (const SteppingHelixPropagator* shp = dynamic_cast<const SteppingHelixPropagator*>(propagator_))
00206      {
00207         SteppingHelixStateInfo state;
00208         try { 
00209            state = shp->propagate(fullTrajectory_[closestPointOnLeft], *plane);
00210         }
00211         catch(cms::Exception &ex){
00212            edm::LogWarning("TrackAssociator") << 
00213                 "Caught exception " << ex.category() << ": " << ex.explainSelf();
00214            edm::LogWarning("TrackAssociator") << "An exception is caught during the track propagation\n"
00215              << state.momentum().x() << ", " << state.momentum().y() << ", " << state.momentum().z();
00216            return TrajectoryStateOnSurface();
00217         }
00218         return state.getStateOnSurface(*plane);
00219      }
00220    else
00221      {
00222         FreeTrajectoryState fts;
00223         fullTrajectory_[closestPointOnLeft].getFreeState(fts);
00224         return propagator_->propagate(fts, *plane);
00225      }
00226 }

bool CachedTrajectory::propagateAll ( const SteppingHelixStateInfo initialState  ) 

propagate through the whole detector, returns true if successful

Definition at line 74 of file CachedTrajectory.cc.

References lat::endl(), Exception, fullTrajectory_, fullTrajectoryFilled_, SteppingHelixStateInfo::isValid(), LogTrace, PV3DBase< T, PVType, FrameType >::mag(), maxRho_, maxZ_, minRho_, minZ_, PV3DBase< T, PVType, FrameType >::perp(), PV3DBase< T, PVType, FrameType >::phi(), SteppingHelixStateInfo::position(), propagateForward(), propagator_, reset_trajectory(), step_, and PV3DBase< T, PVType, FrameType >::z().

Referenced by TrackDetectorAssociator::associate().

00075 {
00076    if ( fullTrajectoryFilled_ ) {
00077       edm::LogWarning("TrackAssociator") << "Reseting all trajectories. Please call reset_trajectory() explicitely to avoid this message";
00078       reset_trajectory();
00079    }
00080         
00081 //   TimerStack timers(TimerStack::Disableable);
00082 
00083    reset_trajectory();
00084    if (propagator_==0) throw cms::Exception("FatalError") << "Track propagator is not defined\n";
00085    SteppingHelixStateInfo currentState(initialState);
00086 
00087    while (currentState.position().perp()<maxRho_ && fabs(currentState.position().z())<maxZ_ ){
00088       LogTrace("TrackAssociator") << "[propagateAll] Propagate outward from (rho, r, z, phi) (" << 
00089         currentState.position().perp() << ", " << currentState.position().mag() << ", " <<
00090         currentState.position().z() << ", " << currentState.position().phi() << ")";
00091       propagateForward(currentState,step_);
00092      if (! currentState.isValid() ) {
00093        LogTrace("TrackAssociator") << "Failed to propagate the track; moving on\n";
00094        break;
00095      }
00096       LogTrace("TrackAssociator") << "\treached (rho, r, z, phi) (" << 
00097         currentState.position().perp() << ", " << currentState.position().mag() << ", " <<
00098         currentState.position().z() << ", " << currentState.position().phi() << ")";
00099      fullTrajectory_.push_back(currentState);
00100    }
00101 
00102 
00103    SteppingHelixStateInfo currentState2(initialState);
00104    SteppingHelixStateInfo previousState;
00105    while (currentState2.position().perp()>minRho_ || fabs(currentState2.position().z())>minZ_) {
00106       previousState=currentState2;
00107       propagateForward(currentState2,-step_);
00108       if (! currentState2.isValid() ) {
00109          LogTrace("TrackAssociator") << "Failed to propagate the track; moving on\n";
00110          break;
00111       }
00112       if(previousState.position().perp()- currentState2.position().perp() < 0) { 
00113          LogTrace("TrackAssociator") << "Error: TrackAssociator has propogated the particle past the point of closest approach to IP" << std::endl;
00114          break;
00115       }
00116       LogTrace("TrackAssociator") << "[propagateAll] Propagated inward from (rho, r, z, phi) (" << 
00117         previousState.position().perp() << ", " << previousState.position().mag() << ", " <<
00118         previousState.position().z() << "," << previousState.position().phi() << ") to (" << 
00119         currentState2.position().perp() << ", " << currentState2.position().mag() << ", " <<
00120         currentState2.position().z() << ", " << currentState2.position().phi() << ")";
00121       fullTrajectory_.push_front(currentState2);
00122    }
00123 
00124 
00125 
00126 
00127    // LogTrace("TrackAssociator") << "fullTrajectory_ has " << fullTrajectory_.size() << " states with (R, z):\n";
00128    // for(uint i=0; i<fullTrajectory_.size(); i++) {
00129    //  LogTrace("TrackAssociator") << "state " << i << ": (" << fullTrajectory_[i].position().perp() << ", "
00130    //    << fullTrajectory_[i].position().z() << ")\n";
00131    // }
00132 
00133 
00134 
00135 
00136 
00137    LogTrace("TrackAssociator") << "Done with the track propagation in the detector. Number of steps: " << fullTrajectory_.size();
00138    fullTrajectoryFilled_ = true;
00139    return ! fullTrajectory_.empty();
00140 }

void CachedTrajectory::propagateForward ( SteppingHelixStateInfo state,
float  distance 
)

Definition at line 28 of file CachedTrajectory.cc.

References cms::Exception::category(), cms::Exception::explainSelf(), TrajectoryStateOnSurface::freeState(), SteppingHelixStateInfo::getFreeState(), SteppingHelixStateInfo::momentum(), SteppingHelixStateInfo::position(), funct::pow(), Propagator::propagate(), propagator_, dttmaxenums::r32, funct::sqrt(), target, Vector3DBase< T, FrameTag >::unit(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by getTrajectory(), and propagateAll().

00029 {
00030    // defined a normal plane wrt the particle trajectory direction
00031    // let's hope that I computed the rotation matrix correctly.
00032    GlobalVector vector(state.momentum().unit());
00033    float r21 = 0;
00034    float r22 = vector.z()/sqrt(1-pow(vector.x(),2));
00035    float r23 = -vector.y()/sqrt(1-pow(vector.x(),2));
00036    float r31 = vector.x();
00037    float r32 = vector.y();
00038    float r33 = vector.z();
00039    float r11 = r22*r33-r23*r32;
00040    float r12 = r23*r31;
00041    float r13 = -r22*r31;
00042    
00043    Surface::RotationType rotation(r11, r12, r13,
00044                                   r21, r22, r23,
00045                                   r31, r32, r33);
00046    Surface* target = new Plane(state.position()+vector*distance, rotation);
00047    if( const SteppingHelixPropagator* shp = dynamic_cast<const SteppingHelixPropagator*>(propagator_) )
00048      {
00049         try {
00050            state = shp->propagate(state, *target);
00051         }
00052         catch(cms::Exception &ex){
00053            edm::LogWarning("TrackAssociator") << 
00054                 "Caught exception " << ex.category() << ": " << ex.explainSelf();
00055            edm::LogWarning("TrackAssociator") << "An exception is caught during the track propagation\n"
00056              << state.momentum().x() << ", " << state.momentum().y() << ", " << state.momentum().z();
00057            state = SteppingHelixStateInfo();
00058         }
00059      }
00060    else
00061      {
00062         FreeTrajectoryState fts;
00063         state.getFreeState( fts );
00064         TrajectoryStateOnSurface stateOnSurface = propagator_->propagate(fts, *target);
00065         state = SteppingHelixStateInfo( *(stateOnSurface.freeState()) );
00066      }
00067    delete target;
00068    // LogTrace("TrackAssociator")
00069    // << state.position().mag() << " , "   << state.position().eta() << " , "
00070    // << state.position().phi();
00071 }

void CachedTrajectory::reset_trajectory (  ) 

Definition at line 365 of file CachedTrajectory.cc.

References ecalTrajectory_, fullTrajectory_, fullTrajectoryFilled_, hcalTrajectory_, hoTrajectory_, wideEcalTrajectory_, wideHcalTrajectory_, and wideHOTrajectory_.

Referenced by TrackDetectorAssociator::associate(), CachedTrajectory(), and propagateAll().

00365                                         { 
00366    fullTrajectory_.clear();
00367    ecalTrajectory_.clear();
00368    hcalTrajectory_.clear();
00369    hoTrajectory_.clear();
00370    wideEcalTrajectory_.clear();   
00371    wideHcalTrajectory_.clear();
00372    wideHOTrajectory_.clear();
00373    fullTrajectoryFilled_ = false;
00374 }

void CachedTrajectory::setMaxDetectorLength ( float  l = 2200.  )  [inline]

Definition at line 89 of file CachedTrajectory.h.

References edm::es::l(), and maxZ_.

Referenced by TrackDetectorAssociator::associate(), and CachedTrajectory().

00089 { maxZ_ = l/2.;}

void CachedTrajectory::setMaxDetectorRadius ( float  r = 800.  )  [inline]

Definition at line 88 of file CachedTrajectory.h.

References maxRho_, and r.

Referenced by TrackDetectorAssociator::associate(), and CachedTrajectory().

00088 { maxRho_ = r;}

void CachedTrajectory::setMaxHOLength ( float  l = 2200.  )  [inline]

Definition at line 91 of file CachedTrajectory.h.

References HOmaxZ_, and edm::es::l().

Referenced by TrackDetectorAssociator::associate().

00091 { HOmaxZ_ = l/2.;}

void CachedTrajectory::setMaxHORadius ( float  r = 800.  )  [inline]

Definition at line 90 of file CachedTrajectory.h.

References HOmaxRho_, and r.

Referenced by TrackDetectorAssociator::associate().

00090 { HOmaxRho_ = r;}

void CachedTrajectory::setMinDetectorLength ( float  l = 0.  )  [inline]

Definition at line 93 of file CachedTrajectory.h.

References edm::es::l(), and minZ_.

Referenced by TrackDetectorAssociator::associate(), and CachedTrajectory().

00093 { minZ_ = l/2.;}

void CachedTrajectory::setMinDetectorRadius ( float  r = 0.  )  [inline]

Definition at line 92 of file CachedTrajectory.h.

References minRho_, and r.

Referenced by TrackDetectorAssociator::associate(), and CachedTrajectory().

00092 { minRho_ = r;}

void CachedTrajectory::setPropagationStep ( float  s = 20.  )  [inline]

Definition at line 95 of file CachedTrajectory.h.

References s, and step_.

Referenced by TrackDetectorAssociator::associate(), and CachedTrajectory().

00095 { step_ = s;}

void CachedTrajectory::setPropagator ( const Propagator ptr  )  [inline]

Definition at line 54 of file CachedTrajectory.h.

References propagator_.

Referenced by TrackDetectorAssociator::setPropagator().

00054 {       propagator_ = ptr; }

void CachedTrajectory::setStateAtIP ( const SteppingHelixStateInfo state  )  [inline]

Definition at line 55 of file CachedTrajectory.h.

References stateAtIP_.

Referenced by TrackDetectorAssociator::associate().

00055 { stateAtIP_ = state; }

static int CachedTrajectory::sign ( float  number  )  [inline, static, protected]

Definition at line 100 of file CachedTrajectory.h.

Referenced by propagate().

00100                                  {
00101       if (number ==0) return 0;
00102       if (number > 0)
00103         return 1;
00104       else
00105         return -1;
00106    }

std::pair< float, float > CachedTrajectory::trajectoryDelta ( TrajectorType  trajectoryType  ) 

calculate trajectory change (Theta,Phi) delta = final - original

Definition at line 228 of file CachedTrajectory.cc.

References delta(), ecalTrajectory_, FullTrajectory, fullTrajectory_, hcalTrajectory_, hoTrajectory_, IpToEcal, IpToHcal, IpToHO, SteppingHelixStateInfo::isValid(), SteppingHelixStateInfo::momentum(), PV3DBase< T, PVType, FrameType >::phi(), HLT_VtxMuL3::result, stateAtIP_, and PV3DBase< T, PVType, FrameType >::theta().

Referenced by TrackDetectorAssociator::fillCaloTowers(), TrackDetectorAssociator::fillEcal(), TrackDetectorAssociator::fillHcal(), TrackDetectorAssociator::fillHO(), and TrackDetectorAssociator::getTAMuonChamberMatches().

00229 {
00230    // MEaning of trajectory change depends on its usage. In most cases we measure 
00231    // change in a trajectory as difference between final track position and initial 
00232    // direction. In some cases such as change of trajectory in the muon detector we 
00233    // might want to compare theta-phi of two points or even find local maximum and
00234    // mimimum. In general it's not essential what defenition of the trajectory change
00235    // is used since we use these numbers only as a rough estimate on how much wider
00236    // we should make the preselection region.
00237    std::pair<float,float> result(0,0);
00238    if ( ! stateAtIP_.isValid() ) { 
00239       edm::LogWarning("TrackAssociator") << "State at IP is not known set. Cannot estimate trajectory change. " <<
00240         "Trajectory change is not taken into account in matching";
00241       return result;
00242    }
00243    switch (trajectoryType) {
00244     case IpToEcal:
00245       if ( ecalTrajectory_.empty() )
00246         edm::LogWarning("TrackAssociator") << "ECAL trajector is empty. Cannot estimate trajectory change. " <<
00247         "Trajectory change is not taken into account in matching";
00248       else return delta( stateAtIP_.momentum().theta(), ecalTrajectory_.front().position().theta(), 
00249                          stateAtIP_.momentum().phi(), ecalTrajectory_.front().position().phi() );
00250       break;
00251     case IpToHcal:
00252       if ( hcalTrajectory_.empty() )
00253         edm::LogWarning("TrackAssociator") << "HCAL trajector is empty. Cannot estimate trajectory change. " <<
00254         "Trajectory change is not taken into account in matching";
00255       else return delta( stateAtIP_.momentum().theta(), hcalTrajectory_.front().position().theta(), 
00256                          stateAtIP_.momentum().phi(),   hcalTrajectory_.front().position().phi() );
00257       break;
00258     case IpToHO:
00259       if ( hoTrajectory_.empty() )
00260         edm::LogWarning("TrackAssociator") << "HO trajector is empty. Cannot estimate trajectory change. " <<
00261         "Trajectory change is not taken into account in matching";
00262       else return delta( stateAtIP_.momentum().theta(), hoTrajectory_.front().position().theta(), 
00263                          stateAtIP_.momentum().phi(),   hoTrajectory_.front().position().phi() );
00264       break;
00265     case FullTrajectory:
00266       if ( fullTrajectory_.empty() )
00267         edm::LogWarning("TrackAssociator") << "Full trajector is empty. Cannot estimate trajectory change. " <<
00268         "Trajectory change is not taken into account in matching";
00269       else  return delta( stateAtIP_.momentum().theta(), fullTrajectory_.back().position().theta(), 
00270                           stateAtIP_.momentum().phi(),   fullTrajectory_.back().position().phi() );
00271       break;
00272     default:
00273       edm::LogWarning("TrackAssociator") << "Unkown or not supported trajector type. Cannot estimate trajectory change. " <<
00274         "Trajectory change is not taken into account in matching";
00275    }
00276    return result;
00277 }


Member Data Documentation

std::vector<SteppingHelixStateInfo> CachedTrajectory::ecalTrajectory_ [protected]

Definition at line 119 of file CachedTrajectory.h.

Referenced by findEcalTrajectory(), getEcalTrajectory(), getStateAtEcal(), reset_trajectory(), and trajectoryDelta().

std::deque<SteppingHelixStateInfo> CachedTrajectory::fullTrajectory_ [protected]

Definition at line 118 of file CachedTrajectory.h.

Referenced by distance(), getInnerState(), getOuterState(), getTrajectory(), propagate(), propagateAll(), reset_trajectory(), and trajectoryDelta().

bool CachedTrajectory::fullTrajectoryFilled_ [protected]

Definition at line 127 of file CachedTrajectory.h.

Referenced by getTrajectory(), propagateAll(), and reset_trajectory().

std::vector<SteppingHelixStateInfo> CachedTrajectory::hcalTrajectory_ [protected]

Definition at line 120 of file CachedTrajectory.h.

Referenced by findHcalTrajectory(), getHcalTrajectory(), getStateAtHcal(), reset_trajectory(), and trajectoryDelta().

float CachedTrajectory::HOmaxRho_ [protected]

Definition at line 133 of file CachedTrajectory.h.

Referenced by setMaxHORadius().

float CachedTrajectory::HOmaxZ_ [protected]

Definition at line 134 of file CachedTrajectory.h.

Referenced by setMaxHOLength().

std::vector<SteppingHelixStateInfo> CachedTrajectory::hoTrajectory_ [protected]

Definition at line 121 of file CachedTrajectory.h.

Referenced by findHOTrajectory(), getHOTrajectory(), getStateAtHO(), reset_trajectory(), and trajectoryDelta().

float CachedTrajectory::maxRho_ [protected]

Definition at line 131 of file CachedTrajectory.h.

Referenced by propagateAll(), and setMaxDetectorRadius().

float CachedTrajectory::maxZ_ [protected]

Definition at line 132 of file CachedTrajectory.h.

Referenced by propagateAll(), and setMaxDetectorLength().

float CachedTrajectory::minRho_ [protected]

Definition at line 135 of file CachedTrajectory.h.

Referenced by propagateAll(), and setMinDetectorRadius().

float CachedTrajectory::minZ_ [protected]

Definition at line 136 of file CachedTrajectory.h.

Referenced by propagateAll(), and setMinDetectorLength().

const Propagator* CachedTrajectory::propagator_ [protected]

Definition at line 129 of file CachedTrajectory.h.

Referenced by propagate(), propagateAll(), propagateForward(), and setPropagator().

SteppingHelixStateInfo CachedTrajectory::stateAtIP_ [protected]

Definition at line 125 of file CachedTrajectory.h.

Referenced by setStateAtIP(), and trajectoryDelta().

float CachedTrajectory::step_ [protected]

Definition at line 137 of file CachedTrajectory.h.

Referenced by getPropagationStep(), propagateAll(), and setPropagationStep().

std::vector<GlobalPoint> CachedTrajectory::wideEcalTrajectory_ [protected]

Definition at line 122 of file CachedTrajectory.h.

Referenced by getWideTrajectory(), and reset_trajectory().

std::vector<GlobalPoint> CachedTrajectory::wideHcalTrajectory_ [protected]

Definition at line 123 of file CachedTrajectory.h.

Referenced by getWideTrajectory(), and reset_trajectory().

std::vector<GlobalPoint> CachedTrajectory::wideHOTrajectory_ [protected]

Definition at line 124 of file CachedTrajectory.h.

Referenced by getWideTrajectory(), and reset_trajectory().


The documentation for this class was generated from the following files:
Generated on Tue Jun 9 18:15:36 2009 for CMSSW by  doxygen 1.5.4