CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
HICMuonUpdator Class Reference

#include <HICMuonUpdator.h>

Public Member Functions

 HICMuonUpdator (double &la1, double &la2, const MagneticField *mf, const cms::HICConst *hh)
 
TrajectoryStateOnSurface update (const Trajectory &mt, const TrajectoryStateOnSurface &, const TrajectoryMeasurement &, const DetLayer *, double &, double &) const
 
TrajectoryStateOnSurface updateBarrel (std::vector< double > &rhit, std::vector< double > &zhit, std::vector< double > &dphihit, std::vector< double > &drhit, std::vector< double > &ehitstrip, std::vector< double > &dehitphi, const TransientTrackingRecHit::ConstRecHitPointer &pRecHit, const TransientTrackingRecHit::ConstRecHitPointer &nRecHit, const TrajectoryStateOnSurface &nTsos, double &, double &, int &) const
 
TrajectoryStateOnSurface updateEndcap (std::vector< double > &rhit, std::vector< double > &zhit, std::vector< double > &dphihit, std::vector< double > &drhit, std::vector< double > &ehitstrip, std::vector< double > &dehitphi, const TransientTrackingRecHit::ConstRecHitPointer &pRecHit, const TransientTrackingRecHit::ConstRecHitPointer &nRecHit, const TrajectoryStateOnSurface &nTsos, double &, double &, int &) const
 
virtual ~HICMuonUpdator ()
 

Private Member Functions

double findPhiInVertex (const FreeTrajectoryState &fts, const double &rc, const GeomDetEnumerators::Location) const
 
bool linefit1 (const std::vector< double > &x, const std::vector< double > &y, const std::vector< double > &err, double &a, double &chi) const
 
bool linefit2 (const std::vector< double > &x, const std::vector< double > &y, const std::vector< double > &err, double &a, double &b, double &chi) const
 

Private Attributes

const MagneticField * field
 
const cms::HICConst * theHICConst
 
double thePhiWin
 
double theZWin
 
double zvert
 

Detailed Description

Definition at line 13 of file HICMuonUpdator.h.

Constructor & Destructor Documentation

HICMuonUpdator::HICMuonUpdator ( double &  la1,
double &  la2,
const MagneticField *  mf,
const cms::HICConst *  hh 
)
inline

Definition at line 16 of file HICMuonUpdator.h.

References field, theHICConst, thePhiWin, theZWin, cms::HICConst::zvert, and zvert.

16  {theHICConst=hh; zvert=hh->zvert;
17  thePhiWin=la1; theZWin=la2; field = mf;};
cms::HICConst::zvert
float zvert
Definition: HICConst.h:25
HICMuonUpdator::field
const MagneticField * field
Definition: HICMuonUpdator.h:49
HICMuonUpdator::theHICConst
const cms::HICConst * theHICConst
Definition: HICMuonUpdator.h:50
virtual HICMuonUpdator::~HICMuonUpdator ( )
inlinevirtual

Definition at line 18 of file HICMuonUpdator.h.

18 {}

Member Function Documentation

double HICMuonUpdator::findPhiInVertex ( const FreeTrajectoryState &  fts,
const double &  rc,
const GeomDetEnumerators::Location  location 
) const
private

Definition at line 315 of file HICMuonUpdator.cc.

References GeomDetEnumerators::barrel, GlobalTrajectoryParameters::charge(), gather_cfg::cout, GlobalTrajectoryParameters::momentum(), FreeTrajectoryState::parameters(), PV3DBase< T, PVType, FrameType >::perp(), PV3DBase< T, PVType, FrameType >::phi(), pi, pi2, GlobalTrajectoryParameters::position(), and PV3DBase< T, PVType, FrameType >::z().

315  {
316 
317  double pi = 4.*atan(1.);
318 // double twopi=8.*atan(1.);
319 
320  double acharge=fts.parameters().charge();
321  double phiclus=fts.parameters().position().phi();
322  double psi;
323  if(location==GeomDetEnumerators::barrel){
324  double xrclus=fts.parameters().position().perp();
325  double xdouble=xrclus/(2.*rc);
326  psi= phiclus+acharge*asin(xdouble);
327  } else {
328  double zclus=fts.parameters().position().z();
329  double pl=fts.parameters().momentum().z();
330  psi=phiclus+acharge*0.006*fabs(zclus)/fabs(pl);
331  }
332  double phic = psi-acharge*pi/2.;
333 #ifdef CORRECT_DEBUG
334  cout<<"Momentum of track="<<fts.parameters().momentum().perp()<<
335  " rad of previous cluster= "<<fts.parameters().position().perp()<<
336  " phi of previous cluster="<<fts.parameters().position().phi()<<endl;
337  cout<<" position of the previous cluster="<<fts.parameters().position()<<endl;
338  cout<<"radius of track="<<rc<<endl;
339  cout<<"acharge="<<acharge<<endl;
340  cout<<"psi="<<psi<<endl;
341  cout<<"phic="<<phic<<" pi="<<pi<<" pi2="<<pi2<<endl;
342 #endif
343 
344  return phic;
345 }
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:71
FreeTrajectoryState::parameters
const GlobalTrajectoryParameters & parameters() const
Definition: FreeTrajectoryState.h:100
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:68
pi2
const double pi2
Definition: Thrust.cc:5
psi
std::map< std::string, int, std::less< std::string > > psi
Definition: CountProcessesAction.h:16
PV3DBase::z
T z() const
Definition: PV3DBase.h:63
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
pi
double pi
Definition: Vx3DHLTAnalyzer.h:63
bool HICMuonUpdator::linefit1 ( const std::vector< double > &  x,
const std::vector< double > &  y,
const std::vector< double > &  err,
double &  a,
double &  chi 
) const
private

Definition at line 270 of file HICMuonUpdator.cc.

References gather_cfg::cout, i, and indexGen::s2.

272  {
273 double s1=0.;
274 double s2=0.;
275 
276 bool fit;
277 fit=false;
278 
279 #ifdef LINEFIT_DEBUG
280 cout<<" linefit1::sizes="<<x.size()<<" "<<y.size()<<" "<<err.size()<<endl;
281 #endif
282 
283 if(x.size() != y.size()) return fit;
284 if(x.size() != err.size()) return fit;
285 
286 for (unsigned int i=0;i<x.size();i++){
287 s1=s1+(x[i]/err[i])*(y[i]/err[i]);
288 s2=s2+(x[i]/err[i])*(x[i]/err[i]);
289 
290 #ifdef LINEFIT_DEBUG
291 cout<<"linefit1="<<x[i]<<" "<<y[i]<<" "<<err[i]<<" "<<s1<<" "<<s2<<endl;
292 #endif
293 
294 }
295 
296 co1=s1/s2;
297 
298 #ifdef LINEFIT_DEBUG
299 cout<<"linefit1::co1,co2="<<co1<<endl;
300 #endif
301 // CHI2
302 chi2=0.;
303 for (unsigned int i=0;i<x.size();i++){
304 chi2=chi2+(y[i]-co1*x[i])*(y[i]-co1*x[i])/(err[i]*err[i]);
305 }
306 chi2=chi2/x.size();
307 #ifdef LINEFIT_DEBUG
308 cout<<"linefit1::chi2="<<chi2<<endl;
309 #endif
310 return fit=true;
311 }
i
int i
Definition: DBlmapReader.cc:9
indexGen.s2
tuple s2
Definition: indexGen.py:106
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
x
Definition: DDAxes.h:10
bool HICMuonUpdator::linefit2 ( const std::vector< double > &  x,
const std::vector< double > &  y,
const std::vector< double > &  err,
double &  a,
double &  b,
double &  chi 
) const
private

Definition at line 170 of file HICMuonUpdator.cc.

References gather_cfg::cout, and i.

172  {
173 double a00=0.;
174 double a01=0.;
175 double a10=0.;
176 double a11=0.;
177 
178 double b0=0.;
179 double b1=0.;
180 
181 bool fit;
182 fit=false;
183 
184 #ifdef LINEFIT_DEBUG
185 cout<<" linefit2::sizes="<<x.size()<<" "<<y.size()<<" "<<err.size()<<endl;
186 cout<<" linefit2::a00 "<<a00<<" "<<a01<<" "<<a11<<" "<<b0<<" "<<b1<<endl;
187 #endif
188 
189 if(x.size() != y.size()) return fit;
190 if(x.size() != err.size()) return fit;
191 
192 
193 for (unsigned int i=0;i<x.size();i++){
194 #ifdef LINEFIT_DEBUG
195 cout<<" line2fit "<<a00<<" "<<x[i]/err[i]<<" "<<x[i]<<" "<<err[i]<<" second try "<<err[i]<<endl;
196 #endif
197 a00=a00+(x[i]/err[i])*(x[i]/err[i]);
198 a01=a01+(x[i]/err[i])/err[i];
199 a10=a01;
200 a11=a11+1./(err[i]*err[i]);
201 b0=b0+(x[i]/err[i])*(y[i]/err[i]);
202 b1=b1+(y[i]/err[i])/err[i];
203 
204 #ifdef LINEFIT_DEBUG
205 cout<<"linefit2="<<x[i]<<" "<<y[i]<<" "<<err[i]<<" "<<a00<<" "<<a01<<" "<<a11<<" "<<b0<<" "<<b1<<endl;
206 #endif
207 
208 }
209 double det=a00*a11-a01*a01;
210 
211 #ifdef LINEFIT_DEBUG
212 cout<<" linefit2::det="<<det<<endl;
213 #endif
214 
215 if(fabs(det)<0.00000001) return fit;
216 co1=(b0*a11-b1*a01)/det;
217 co2=(b1*a00-b0*a10)/det;
218 
219 // check if it is 90 degree track
220 #ifdef LINEFIT_DEBUG
221 cout<<" linefit2::Previous element= "<<y[x.size()-3]<<" "<<x[x.size()-3]<<endl;
222 #endif
223 
224 if(y[x.size()-2]<14.) {
225 
226 #ifdef LINEFIT_DEBUG
227 cout<<" check 90 degree track "<<endl;
228 #endif
229 
230 if(fabs(x[x.size()-2]-x[x.size()-1])<0.1){
231 
232 #ifdef LINEFIT_DEBUG
233 cout<<" Redetermine line - 90 degree "<<endl;
234 #endif
235 
236 if(fabs(x[x.size()-2]-x[x.size()-1])>0.0001){
237 co1=(y[x.size()-2]-y[x.size()-1])/(x[x.size()-2]-x[x.size()-1]);
238 co2=y[x.size()-2]-co1*x[x.size()-2];
239 }
240  else
241  {
242 co1=10000.;
243 co2=-zvert*10000.;
244  }
245 }
246 }
247 
248 #ifdef LINEFIT_DEBUG
249 cout<<"linefit2::co1,co2="<<co1<<" "<<co2<<" size "<<x.size()<<endl;
250 #endif
251 // CHI2
252 chi2=0.;
253 for (unsigned int i=0;i<x.size();i++){
254 double zdet=(y[i]-co2)/co1;
255 chi2=chi2+(x[i]-zdet)*(x[i]-zdet)/(err[i]*err[i]);
256 
257 #ifdef LINEFIT_DEBUG
258 cout<<"linefit2::chi2="<<chi2<<" err="<<err[i]<<" x[i]="<<x[i]<<" teor="<<zdet<<endl;
259 #endif
260 
261 }
262 // chi2 on degree of freedom
263 chi2=chi2/x.size();
264 #ifdef LINEFIT_DEBUG
265 cout<<" linefit2::chi2= "<<chi2<<endl;
266 #endif
267 
268 return fit=true;
269 }
i
int i
Definition: DBlmapReader.cc:9
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
x
Definition: DDAxes.h:10
TrajectoryStateOnSurface HICMuonUpdator::update ( const Trajectory &  mt,
const TrajectoryStateOnSurface &  nTsos,
const TrajectoryMeasurement &  ntm,
const DetLayer *  layer,
double &  chirz,
double &  chirf 
) const

Definition at line 23 of file HICMuonUpdator.cc.

References GeomDetEnumerators::barrel, FreeTrajectoryState::charge(), gather_cfg::cout, TrajectoryStateOnSurface::freeTrajectoryState(), TrajectoryStateOnSurface::isValid(), TrajectoryMeasurement::layer(), DetLayer::location(), Trajectory::measurements(), pi, TrajectoryMeasurement::recHit(), mathSSE::sqrt(), funct::tan(), PV3DBase< T, PVType, FrameType >::theta(), and theta().

Referenced by progressbar.ProgressBar::__next__(), relval_steps.Matrix::__setitem__(), relval_steps.Steps::__setitem__(), Vispa.Gui.VispaWidget.VispaWidget::autosize(), Vispa.Views.LineDecayView.LineDecayContainer::createObject(), Vispa.Views.LineDecayView.LineDecayContainer::deselectAllObjects(), Vispa.Gui.VispaWidgetOwner.VispaWidgetOwner::deselectAllWidgets(), Vispa.Gui.VispaWidget.VispaWidget::enableAutosizing(), progressbar.ProgressBar::finish(), Vispa.Gui.MenuWidget.MenuWidget::leaveEvent(), Vispa.Gui.VispaWidgetOwner.VispaWidgetOwner::mouseMoveEvent(), Vispa.Gui.MenuWidget.MenuWidget::mouseMoveEvent(), Vispa.Views.LineDecayView.LineDecayContainer::mouseMoveEvent(), Vispa.Gui.VispaWidgetOwner.VispaWidgetOwner::mouseReleaseEvent(), Vispa.Views.LineDecayView.LineDecayContainer::objectMoved(), relval_steps.Steps::overwrite(), Vispa.Views.LineDecayView.LineDecayContainer::removeObject(), Vispa.Gui.ConnectableWidget.ConnectableWidget::removePorts(), Vispa.Gui.FindDialog.FindDialog::reset(), Vispa.Gui.PortConnection.PointToPointConnection::select(), Vispa.Gui.VispaWidget.VispaWidget::select(), Vispa.Views.LineDecayView.LineDecayContainer::select(), Vispa.Gui.VispaWidget.VispaWidget::setText(), Vispa.Gui.VispaWidget.VispaWidget::setTitle(), Vispa.Gui.ZoomableWidget.ZoomableWidget::setZoom(), Vispa.Views.LineDecayView.LineDecayContainer::setZoom(), Vispa.Gui.PortConnection.PointToPointConnection::updateConnection(), and HICTrajectoryBuilder::updateTrajectory().

27  {
28  double pi = 4.*atan(1.);
29  double twopi=8.*atan(1.);
30  TrajectoryStateOnSurface badtsos;
31  if(!nTsos.isValid()) {
32  std::cout<<" HICMuonUpdator::update:: can not start::initial tsos is not valid " <<std::endl;
33  return badtsos;
34  }
35 // trajectory type
36 
37  vector<TrajectoryMeasurement> MTM=mt.measurements();
38 #ifdef DEBUG
39  std::cout<<" HICMuonUpdator::update::MTM size "<<MTM.size()<<" vertex "<<zvert<<std::endl;
40  std::cout<<" HICMuonUpdator::update::charge from trajectory"<<(MTM.back()).updatedState().freeTrajectoryState()->parameters().charge()<<std::endl;
41  std::cout<<" HICMuonUpdator::update::charge predicted tsos"<<nTsos.freeTrajectoryState()->charge()<<std::endl;
42  std::cout<<" HICMuonUpdator::update::momentum "<<(MTM.back()).updatedState().freeTrajectoryState()->parameters().momentum()<<std::endl;
43 #endif
44  vector<double> phihit,rhit,zhit,dphihit,drhit,dzhit,dzhitl,ehitphi,dehitphi,ehitstrip;
45 
46  double rvert=0.;
47  double ezvert=0.014;
48 
49  const TransientTrackingRecHit::ConstRecHitPointer pRecHit=(MTM.back()).recHit();
50  const TransientTrackingRecHit::ConstRecHitPointer nRecHit=ntm.recHit();
51 
52 // double acharge=(MTM.back()).updatedState().freeTrajectoryState()->parameters().charge();
53 // double acharge=nTsos.freeTrajectoryState()->charge();
54  GlobalVector pold=(MTM.back()).updatedState().freeTrajectoryState()->parameters().momentum();
55  double theta=pold.theta();
56 
57  for(vector<TrajectoryMeasurement>::const_iterator ihit=MTM.begin();ihit!=MTM.end();ihit++){
58 
59 // FreeTrajectoryState* ftshit = (*ihit).updatedState().freeTrajectoryState();
60  phihit.push_back((*ihit).recHit()->globalPosition().phi());
61  rhit.push_back((*ihit).recHit()->globalPosition().perp());
62  zhit.push_back((*ihit).recHit()->globalPosition().z());
63  GlobalPoint realhit = (*ihit).recHit()->globalPosition();
64 
65  double phierror=sqrt((*ihit).recHit()->globalPositionError().phierr(realhit));
66 
67  if(fabs(phierror)<0.0000001) {
68  phierror=0.00008;
69  }
70  ehitphi.push_back(phierror);
71 
72 #ifdef UPDATOR_BARREL_DEBUG
73  cout<<" Errors "<<phierror<<" "<<(*ihit).recHit()->globalPositionError().rerr(realhit)<<" "<<tan(theta)<<endl;
74 #endif
75 
76 
77  if((*ihit).layer()->location()==GeomDetEnumerators::barrel){
78  ehitstrip.push_back(sqrt((*ihit).recHit()->globalPositionError().czz()));
79  } else{
80  ehitstrip.push_back(sqrt((*ihit).recHit()->globalPositionError().rerr(realhit)/fabs(tan(theta))));
81  }
82 
83  }
84 
85  phihit.push_back(nRecHit->globalPosition().phi());
86  rhit.push_back(nRecHit->globalPosition().perp());
87  zhit.push_back(nRecHit->globalPosition().z());
88  ehitphi.push_back(sqrt(nRecHit->globalPositionError().phierr(nRecHit->globalPosition())));
89 
90 
91 
92 
93 
94  if(ntm.layer()->location()==GeomDetEnumerators::barrel){
95  ehitstrip.push_back(sqrt(nRecHit->globalPositionError().czz()));
96  } else {
97  ehitstrip.push_back(sqrt(nRecHit->globalPositionError().rerr(nRecHit->globalPosition()))/fabs(tan(theta)));
98  }
99 
100 // add vertex
101 
102  rhit.push_back(rvert);
103  zhit.push_back(zvert);
104  ehitstrip.push_back(ezvert);
105 
106  for(vector<double>::const_iterator iphi=phihit.begin();iphi!=phihit.end()-1;iphi++){
107  double dpnew=fabs(*iphi-*(iphi+1));
108  if(dpnew>pi) dpnew=twopi-dpnew;
109 
110 #ifdef UPDATOR_BARREL_DEBUG
111  cout<<" dphi=dpnew="<<dpnew<<" "<<*iphi<<" "<<*(iphi+1)<<endl;
112 #endif
113 
114  dphihit.push_back(dpnew);
115  }
116 
117  for(vector<double>::const_iterator ir=rhit.begin();ir!=rhit.end()-2;ir++){
118  double dpnew=fabs(*ir-*(ir+1));
119 
120 #ifdef UPDATOR_BARREL_DEBUG
121  cout<<" dr=dpnew="<<dpnew<<" "<<*ir<<" "<<*(ir+1)<<endl;
122 #endif
123 
124  drhit.push_back(dpnew);
125  }
126  for(vector<double>::const_iterator iz=zhit.begin();iz!=zhit.end()-2;iz++){
127  double dpnew=*iz-*(iz+1);
128 
129 #ifdef UPDATOR_BARREL_DEBUG
130  cout<<" dZ=dpnew="<<dpnew<<" "<<*iz<<" "<<*(iz+1)<<endl;
131 #endif
132 
133  dzhit.push_back(fabs(dpnew));
134  dzhitl.push_back(dpnew);
135  }
136 
137  dzhitl.push_back(*(zhit.end()-1)-zvert);
138 
139  for(vector<double>::const_iterator ie=ehitphi.begin();ie!=ehitphi.end()-1;ie++){
140  double dpnew=(*ie)*1.14;
141  dehitphi.push_back(dpnew);
142  }
143 //
144 //=================fit in rf and rz planes separately
145 //
146 int tType = 1;
147 if ( (*(MTM.begin())).layer()->location()==GeomDetEnumerators::barrel){
148 // std::cout<<" Update barrel "<<std::endl;
149  TrajectoryStateOnSurface tsos=updateBarrel(rhit, zhit, dphihit, drhit, ehitstrip, dehitphi, pRecHit, nRecHit,
150  nTsos, chirz, chirf, tType);
151 
152  if(!tsos.isValid()) {
153 #ifdef DEBUG
154  std::cout<<" Trajectory is not valid "<<std::endl;
155 #endif
156  return badtsos;
157  }
158 
159  return tsos;
160 
161 } else{
162 // std::cout<<" Update endcap "<<std::endl;
163  TrajectoryStateOnSurface tsos=updateEndcap(rhit, zhit, dphihit, dzhit, ehitstrip, dehitphi, pRecHit, nRecHit,
164  nTsos, chirz, chirf, tType);
165  return tsos;
166 
167 }
168 }
DetLayer::location
virtual Location location() const =0
Which part of the detector (barrel, endcap)
TrajectoryStateOnSurface::freeTrajectoryState
FreeTrajectoryState * freeTrajectoryState(bool withErrors=true) const
Definition: TrajectoryStateOnSurface.h:134
theta
Geom::Theta< T > theta() const
Definition: Basic3DVectorLD.h:165
TrajectoryMeasurement::recHit
ConstRecHitPointer recHit() const
Definition: TrajectoryMeasurement.h:191
FreeTrajectoryState::charge
TrackCharge charge() const
Definition: FreeTrajectoryState.h:81
PV3DBase::theta
Geom::Theta< T > theta() const
Definition: PV3DBase.h:74
Trajectory::measurements
DataContainer const & measurements() const
Definition: Trajectory.h:203
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:46
funct::tan
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
TrajectoryMeasurement::layer
const DetLayer * layer() const
Definition: TrajectoryMeasurement.h:197
HICMuonUpdator::updateBarrel
TrajectoryStateOnSurface updateBarrel(std::vector< double > &rhit, std::vector< double > &zhit, std::vector< double > &dphihit, std::vector< double > &drhit, std::vector< double > &ehitstrip, std::vector< double > &dehitphi, const TransientTrackingRecHit::ConstRecHitPointer &pRecHit, const TransientTrackingRecHit::ConstRecHitPointer &nRecHit, const TrajectoryStateOnSurface &nTsos, double &, double &, int &) const
Definition: HICMuonUpdator.cc:347
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
pi
double pi
Definition: Vx3DHLTAnalyzer.h:63
HICMuonUpdator::updateEndcap
TrajectoryStateOnSurface updateEndcap(std::vector< double > &rhit, std::vector< double > &zhit, std::vector< double > &dphihit, std::vector< double > &drhit, std::vector< double > &ehitstrip, std::vector< double > &dehitphi, const TransientTrackingRecHit::ConstRecHitPointer &pRecHit, const TransientTrackingRecHit::ConstRecHitPointer &nRecHit, const TrajectoryStateOnSurface &nTsos, double &, double &, int &) const
Definition: HICMuonUpdator.cc:445
TrajectoryStateOnSurface HICMuonUpdator::updateBarrel ( std::vector< double > &  rhit,
std::vector< double > &  zhit,
std::vector< double > &  dphihit,
std::vector< double > &  drhit,
std::vector< double > &  ehitstrip,
std::vector< double > &  dehitphi,
const TransientTrackingRecHit::ConstRecHitPointer &  pRecHit,
const TransientTrackingRecHit::ConstRecHitPointer &  nRecHit,
const TrajectoryStateOnSurface &  nTsos,
double &  chirz,
double &  chirf,
int &  tType 
) const

Definition at line 347 of file HICMuonUpdator.cc.

References GlobalTrajectoryParameters::charge(), funct::cos(), gather_cfg::cout, TrajectoryStateOnSurface::freeTrajectoryState(), m, FreeTrajectoryState::parameters(), PV3DBase< T, PVType, FrameType >::perp(), PV3DBase< T, PVType, FrameType >::phi(), pi, funct::sin(), and TrajectoryStateOnSurface::surface().

353  {
354 
355 // fit in (dphi dr), (dphi-dz)
356  TrajectoryStateOnSurface badtsos;
357  double pi = 4.*atan(1.);
358  double twopi=8.*atan(1.);
359 
360 // cout<<" Update barrel begin "<<endl;
361 
362  double ch1,dphi,dr,ptnew;
363  double co1,co2;
364  bool fitrf,fitrz;
365 
366 // fit in (ZR)-coordinate
367 
368  fitrz=this->linefit2(rhit,zhit,ehitstrip,co1,co2,chirz);
369 
370 #ifdef UPDATOR_BARREL_DEBUG
371  cout<<"UPDATE::BARREL::line fit rz= "<<fitrz<<" chirz="<<chirz<<endl;
372  cout<<" co1="<<co1<<" co2="<<co2<<endl;
373 #endif
374 
375  if(!fitrz) {
376 #ifdef DEBUG
377  cout<<"UPDATE::BARREL::line fit failed rz= "<<fitrz<<" chirz="<<chirz<<endl;
378 #endif
379 
380  return badtsos;
381  }
382  if(dphihit.size()>1){
383  fitrf=this->linefit1(dphihit,drhit,dehitphi,ch1,chirf);
384 
385 #ifdef UPDATOR_BARREL_DEBUG
386  cout<<"UPDATE::BARREL::line fit dphi= "<<fitrf<<" chirf="<<chirf<<endl;
387  cout<<" ch1="<<ch1<<endl;
388 #endif
389 
390  if(!fitrf) {
391 #ifdef DEBUG
392  cout<<"UPDATE::BARREL::line fit failed dphi= "<<fitrf<<" chirz="<<chirf<<endl;
393 #endif
394  return badtsos;
395  }
396  }else{
397 
398  chirf = 0.;
399  dphi=fabs(dphihit.back());
400  dr=fabs(drhit.back());
401  if(dphi > pi) dphi = twopi-dphi;
402  ch1=dphi/dr;
403 
404 #ifdef UPDATOR_BARREL_DEBUG
405  cout<<"UPDATE::BARREL::line calc dphi= "<<dphi<<" dr="<<dr<<" chirf="<<chirf<<endl;
406  cout<<" ch1="<<ch1<<endl;
407 #endif
408 
409  }
410 
411 // Updating trajectory
412  ptnew=0.006/ch1;
413  GlobalPoint xrhit = nRecHit->globalPosition();
414  TrackCharge aCharge = nTsos.freeTrajectoryState()->parameters().charge();
415  double phiclus=xrhit.phi();
416  double xrclus=xrhit.perp();
417  //double xzclus=xrhit.z();
418  double rc=100.*ptnew/(0.3*4.);
419  double xdouble=xrclus/(2.*rc);
420  double psi=phiclus-aCharge*asin(xdouble);
421  double pznew=ptnew/co1;
422  double znew=(xrclus-co2)/co1;
423  double delphinew=fabs(0.006*drhit.back()/ptnew);
424  double phinew=pRecHit->globalPosition().phi()+aCharge*delphinew;
425  GlobalVector pnew(ptnew*cos(psi),ptnew*sin(psi),pznew);
426  GlobalPoint xnew(xrclus*cos(phinew),xrclus*sin(phinew),znew);
427  AlgebraicSymMatrix55 m;
428  m(0,0)=ptnew; m(1,1)=thePhiWin,
429  m(2,2)=theZWin,
430  m(3,3)=thePhiWin,
431  m(4,4)=theZWin;
432 
433 #ifdef UPDATOR_BARREL_DEBUG
434  cout<< "MuUpdator::xnew=" << xnew<<endl;
435 #endif
436 
437  TrajectoryStateOnSurface tsos(
438  GlobalTrajectoryParameters(xnew, pnew, aCharge, field),
439  CurvilinearTrajectoryError(m), nTsos.surface());
440 
441 // cout<<" Update barrel end "<<xnew<<endl;
442 
443  return tsos;
444 }
HICMuonUpdator::linefit1
bool linefit1(const std::vector< double > &x, const std::vector< double > &y, const std::vector< double > &err, double &a, double &chi) const
Definition: HICMuonUpdator.cc:270
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:71
FreeTrajectoryState::parameters
const GlobalTrajectoryParameters & parameters() const
Definition: FreeTrajectoryState.h:100
TrajectoryStateOnSurface::freeTrajectoryState
FreeTrajectoryState * freeTrajectoryState(bool withErrors=true) const
Definition: TrajectoryStateOnSurface.h:134
HICMuonUpdator::linefit2
bool linefit2(const std::vector< double > &x, const std::vector< double > &y, const std::vector< double > &err, double &a, double &b, double &chi) const
Definition: HICMuonUpdator.cc:170
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:68
AlgebraicSymMatrix55
ROOT::Math::SMatrix< double, 5, 5, ROOT::Math::MatRepSym< double, 5 > > AlgebraicSymMatrix55
Definition: AlgebraicROOTObjects.h:21
HICMuonUpdator::field
const MagneticField * field
Definition: HICMuonUpdator.h:49
psi
std::map< std::string, int, std::less< std::string > > psi
Definition: CountProcessesAction.h:16
TrackCharge
int TrackCharge
Definition: TrackCharge.h:4
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
TrajectoryStateOnSurface::surface
const Surface & surface() const
Definition: TrajectoryStateOnSurface.h:182
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
pi
double pi
Definition: Vx3DHLTAnalyzer.h:63
TrajectoryStateOnSurface HICMuonUpdator::updateEndcap ( std::vector< double > &  rhit,
std::vector< double > &  zhit,
std::vector< double > &  dphihit,
std::vector< double > &  drhit,
std::vector< double > &  ehitstrip,
std::vector< double > &  dehitphi,
const TransientTrackingRecHit::ConstRecHitPointer &  pRecHit,
const TransientTrackingRecHit::ConstRecHitPointer &  nRecHit,
const TrajectoryStateOnSurface &  nTsos,
double &  chirz,
double &  chirf,
int &  tType 
) const

Definition at line 445 of file HICMuonUpdator.cc.

References FreeTrajectoryState::charge(), funct::cos(), gather_cfg::cout, TrajectoryStateOnSurface::freeTrajectoryState(), m, PV3DBase< T, PVType, FrameType >::phi(), pi, funct::sin(), TrajectoryStateOnSurface::surface(), and PV3DBase< T, PVType, FrameType >::z().

451  {
452 
453 // fit in (dphi dr), (dphi-dz)
454  double pi = 4.*atan(1.);
455  double twopi=8.*atan(1.);
456 
457 
458  TrajectoryStateOnSurface badtsos;
459 
460 // cout<<" Update endcap begin "<<endl;
461 
462  double ch1,dphi,dr;
463  double co1,co2;
464  bool fitrf,fitrz;
465 
466 #ifdef UPDATOR_ENDCAP_DEBUG
467  cout<<"updateEndcap switched on"<<endl;
468 #endif
469 
470 // fit in (ZR)-coordinate
471 
472  fitrz=this->linefit2(rhit,zhit,ehitstrip,co1,co2,chirz);
473 
474 #ifdef UPDATOR_ENDCAP_DEBUG
475  cout<<"line fit rz= "<<fitrz<<" chirz="<<chirz<<endl;
476  cout<<" co1="<<co1<<" co2="<<co2<<endl;
477 #endif
478 
479  if(!fitrz) {
480 #ifdef DEBUG
481  cout<<"UPDATE::ENDCAP::line fit failed rz= "<<fitrz<<" chirz="<<chirz<<endl;
482 #endif
483  return badtsos;
484  }
485  if(dphihit.size()>1){
486  fitrf=this->linefit1(dphihit,drhit,dehitphi,ch1,chirf);
487  if(zhit.front()<0.) ch1=-1.*ch1;
488 #ifdef UPDATOR_ENDCAP_DEBUG
489  cout<<"MuUpdate::barrel::line fit dphi= "<<fitrf<<" chirf="<<chirf<<endl;
490  cout<<" ch1="<<ch1<<endl;
491 #endif
492  if(!fitrf) {
493 #ifdef DEBUG
494  cout<<"UPDATE::ENDCAP::line fit failed dphi= "<<fitrf<<" chirz="<<chirf<<endl;
495 #endif
496 
497  return badtsos;
498  }
499  }else{
500  dphi=fabs(dphihit.back());
501  dr=fabs(drhit.back());
502  if(dphi > pi) dphi = twopi-dphi;
503  ch1=dphi/dr;
504  if(zhit.front()<0.) ch1=-1.*ch1;
505  chirf = 0.;
506 #ifdef UPDATOR_ENDCAP_DEBUG
507  cout<< "MuUpdate::barrel::linedphi=" << dphi <<" dz= "<< dr << " ch1= " <<ch1<<" pznew= "<<0.006/ch1<<endl;
508 #endif
509  }
510 
511 // Updating trajectory
512  double pznew=0.006/ch1;
513 // cout<<" point 1 "<<endl;
514  GlobalPoint xrhit = nRecHit->globalPosition();
515 // cout<<" point 2 "<<endl;
516  TrackCharge aCharge = nTsos.freeTrajectoryState()->charge();
517 // cout<<" point 3 "<<endl;
518  double phiclus=xrhit.phi();
519 // cout<<" point 4 "<<endl;
520  double xzclus=xrhit.z();
521 // cout<<" point 5 "<<endl;
522  double psi=phiclus-aCharge*0.006*fabs(xzclus-zvert)/fabs(pznew);
523 // cout<<" point 6 "<<endl;
524  double ptnew=pznew*co1;
525 // cout<<" point 7 "<<endl;
526  double xrclus=co1*xzclus+co2;
527 // cout<<" point 8 "<<endl;
528  double delphinew=fabs(0.006*drhit.back()/pznew);
529 // cout<<" point 9 "<<endl;
530  double phinew=pRecHit->globalPosition().phi()+aCharge*delphinew;
531 // cout<<" point 10 "<<endl;
532  GlobalVector pnew(ptnew*cos(psi),ptnew*sin(psi),pznew);
533 // cout<<" point 11 "<<endl;
534  GlobalPoint xnew(xrclus*cos(phinew),xrclus*sin(phinew),xzclus);
535 // OK changes. Start each time from the RealHit cluster.
536 //
537 // GlobalPoint xnew(xrhit.x(),xrhit.y(),xzclus);
538 
539 #ifdef UPDATOR_ENDCAP_DEBUG
540  cout<< "MuUpdator::xnew=" << xnew<<endl;
541 #endif
542 
543 
544  AlgebraicSymMatrix55 m;
545  m(0,0)=pznew; m(1,1)=thePhiWin,
546  m(2,2)=theZWin,
547  m(3,3)=thePhiWin,
548  m(4,4)=theZWin;
549 
550  TrajectoryStateOnSurface tsos(
551  GlobalTrajectoryParameters(xnew, pnew, aCharge, field),
552  CurvilinearTrajectoryError(m), nTsos.surface());
553 
554 // cout<< "Update endcap end "<<endl;
555  return tsos;
556 }
HICMuonUpdator::linefit1
bool linefit1(const std::vector< double > &x, const std::vector< double > &y, const std::vector< double > &err, double &a, double &chi) const
Definition: HICMuonUpdator.cc:270
TrajectoryStateOnSurface::freeTrajectoryState
FreeTrajectoryState * freeTrajectoryState(bool withErrors=true) const
Definition: TrajectoryStateOnSurface.h:134
HICMuonUpdator::linefit2
bool linefit2(const std::vector< double > &x, const std::vector< double > &y, const std::vector< double > &err, double &a, double &b, double &chi) const
Definition: HICMuonUpdator.cc:170
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:68
AlgebraicSymMatrix55
ROOT::Math::SMatrix< double, 5, 5, ROOT::Math::MatRepSym< double, 5 > > AlgebraicSymMatrix55
Definition: AlgebraicROOTObjects.h:21
HICMuonUpdator::field
const MagneticField * field
Definition: HICMuonUpdator.h:49
FreeTrajectoryState::charge
TrackCharge charge() const
Definition: FreeTrajectoryState.h:81
psi
std::map< std::string, int, std::less< std::string > > psi
Definition: CountProcessesAction.h:16
TrackCharge
int TrackCharge
Definition: TrackCharge.h:4
PV3DBase::z
T z() const
Definition: PV3DBase.h:63
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
TrajectoryStateOnSurface::surface
const Surface & surface() const
Definition: TrajectoryStateOnSurface.h:182
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
pi
double pi
Definition: Vx3DHLTAnalyzer.h:63

Member Data Documentation

const MagneticField* HICMuonUpdator::field
private

Definition at line 49 of file HICMuonUpdator.h.

Referenced by HICMuonUpdator().

const cms::HICConst* HICMuonUpdator::theHICConst
private

Definition at line 50 of file HICMuonUpdator.h.

Referenced by HICMuonUpdator().

double HICMuonUpdator::thePhiWin
private

Definition at line 47 of file HICMuonUpdator.h.

Referenced by HICMuonUpdator().

double HICMuonUpdator::theZWin
private

Definition at line 48 of file HICMuonUpdator.h.

Referenced by HICMuonUpdator().

double HICMuonUpdator::zvert
private

Definition at line 46 of file HICMuonUpdator.h.

Referenced by HICMuonUpdator().

Generated for CMSSW Reference Manual by doxygen 1.8.5