Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages

MdcxAddHits Class Reference

#include <MdcxAddHits.h>

List of all members.

Public Member Functions

const HepAList< MdcxHit > & GetAssociates (int i=0)
const HepAList< MdcxHit > & GetAssociates (int i=0)
 MdcxAddHits (HepAList< MdcxHel > &f, const HepAList< MdcxHit > &h, float addit=1.5)
 MdcxAddHits (HepAList< MdcxFittedHel > &f, const HepAList< MdcxHit > &h, float addit=1.5)
 MdcxAddHits (HepAList< MdcxHel > &f, const HepAList< MdcxHit > &h, float addit=1.5)
 MdcxAddHits (HepAList< MdcxFittedHel > &f, const HepAList< MdcxHit > &h, float addit=1.5)
virtual ~MdcxAddHits ()
virtual ~MdcxAddHits ()

Protected Attributes

float addcut
HepAList< MdcxHit > hhhh
HepAList< MdcxHit > hhhh
int kkk
int kkl
HepAList< MdcxHit > listoasses
HepAList< MdcxHit > listoasses
int nadded
int ncalc
int nhits
int ntracks
double sumpull
int thot
HepAList< MdcxHel > tttt
HepAList< MdcxHel > tttt
int tuot

Constructor & Destructor Documentation

MdcxAddHits::MdcxAddHits (  HepAList< MdcxFittedHel > &  f,
const HepAList< MdcxHit > &  h,
float  addit = 1.5
) 
  

00041                                                            :
00042   ncalc(0), nadded(0), sumpull(0.0), thot(0), tuot(0)
00043 {
00044   ntracks = trkl.length();
00045   nhits   = hitl.length();
00046   addcut  = addit;  
00047   // cout << "MdcxAddHits called with " << ntracks << " tracks, " 
00048   //      << nhits << " hits, addit = " << addit << endl;
00049   if ( (ntracks < 1) || (nhits < 1) ) return;
00050   kkk = 0;
00051   kkl = 0;
00052   while(hitl[kkl]) {
00053     hitl[kkl]->SetUsedOnHel(0);
00054     hhhh.append(hitl[kkl++]);
00055   }
00056   while(trkl[kkk]) {
00057     MdcxHel* thel = trkl[kkk];
00058     tttt.append(thel);
00059     const HepAList<MdcxHit>& dchitlist = trkl[kkk++]->XHitList();
00060     thot += dchitlist.length();
00061     kkl = 0;
00062     while (dchitlist[kkl]) dchitlist[kkl++]->SetUsedOnHel(1);
00063   }
00064 }

MdcxAddHits::MdcxAddHits (  HepAList< MdcxHel > &  f,
const HepAList< MdcxHit > &  h,
float  addit = 1.5
) 
  

00066                                                                                            :
00067   ncalc(0), nadded(0), sumpull(0.0)
00068 {
00069   //
00070   //  c-tor designed to work with MdcxHitAdder; all MdcxHit's in hitl are
00071   //  assumed to be fresh and fair game.
00072   //
00073   ntracks = trkl.length();
00074   nhits   = hitl.length();
00075   addcut  = addit;  
00076   // cout << "MdcxAddHits called with " << ntracks << " tracks, " 
00077   //      << nhits << " hits, addit = " << addit << endl;
00078   if ( (ntracks < 1) || (nhits < 1) ) return;
00079   hhhh = hitl;
00080   tttt = trkl;
00081 }

MdcxAddHits::~MdcxAddHits (   )  [virtual]
  

00083                           { 
00084 }

MdcxAddHits::MdcxAddHits (  HepAList< MdcxFittedHel > &  f,
const HepAList< MdcxHit > &  h,
float  addit = 1.5
) 
 

MdcxAddHits::MdcxAddHits (  HepAList< MdcxHel > &  f,
const HepAList< MdcxHit > &  h,
float  addit = 1.5
) 
 

virtual MdcxAddHits::~MdcxAddHits (   )  [virtual]
 

Member Function Documentation

const HepAList<MdcxHit>& MdcxAddHits::GetAssociates (  int  i = 0  ) 
 

const HepAList< MdcxHit > & MdcxAddHits::GetAssociates (  int  whichtrack = 0  ) 
 

zoujh

FIXME epsilon ? yzhang

FIXME 

00086                                                                   {
00087   int debug = MdcxParameters::debug;
00088 
00089   listoasses.removeAll();  
00090   if ((whichtrack >= ntracks) || (whichtrack < 0)) {
00091     cout << "asked for associates of track " << whichtrack
00092       << " while there are only " << ntracks << " tracks in list " << endl;
00093     return listoasses;
00094   }
00095 
00096   double m_2pi = 2.0*M_PI;
00097   MdcxHel* temphel = tttt[whichtrack];
00098 
00099   if(5 == debug) temphel->print();
00100 
00101   double lmax = temphel->Lmax(); 
00102 
00103   //  calc. phim, phip
00104   double gvin = MdcxParameters::firstMdcAxialRadius;//yzhang 2010-05-05 
00105   double gvout = MdcxParameters::maxMdcRadius;
00106   double phiex = 0.38;  
00107   double phim = -M_PI;
00108   double phip = M_PI;
00109   double rmin = fabs(temphel->D0());
00110   double rmax, pc, rc, rhel;
00111   if (temphel->Ominfl()) {
00112     rmax = fabs(temphel->D0() + 2.0/temphel->Omega());
00113     double xc = temphel->Xc();
00114     double yc = temphel->Yc();
00115     pc = atan2(yc,xc);  
00116     rc = fabs(temphel->D0() + 1.0/temphel->Omega());
00117     rhel = fabs(1.0/temphel->Omega());
00118   } else {
00119     rmax = 1000.;
00120     pc = temphel->Phi0();
00121     rc = 1000.;
00122     rhel = 1000.;
00123   }
00124   if(5 == debug) {
00125     std::cout<<__FILE__<<" " << " rmin >?"<<rmin << " gvin "<<gvin
00126       << " rmax >?"<<rmax << " gvout "<<gvout <<  std::endl;
00127   }
00128 
00129   if ((rmin<gvin) && (rmax>gvout)) {
00130     // case A (exiter)
00131     if(5 == debug) std::cout<<__FILE__<<" case A (exiter) "<<std::endl;
00132 
00133     double len = 0;
00134     double lstep = m_2pi*rhel/16.;
00135     if (lstep>10.) lstep = 10.;
00136     double xl = temphel->Xh(len);
00137     double yl = temphel->Yh(len);
00138     double phil = atan2(yl, xl);
00139     double rl   = sqrt(xl*xl + yl*yl);
00140     int nstep = 0;
00141     double philast = phil;
00142     int isin   = 0;
00143     int notout = 1;
00144     while ((notout) && (nstep<20)) {
00145       len += lstep;
00146       nstep++;
00147       xl = temphel->Xh(len);
00148       yl = temphel->Yh(len);
00149       phil = atan2(yl, xl);
00150       rl = sqrt(xl*xl + yl*yl);
00151       if ((rl<gvin) && (!isin)) {
00152         philast = phil;
00153       } else if ((rl>gvin) && (!isin)) {
00154         isin = 1; phim = philast; phip = philast;
00155       }
00156       if (isin) {
00157         double deltap = phil - philast;
00158         if (deltap >  M_PI) phil -= m_2pi;
00159         if (deltap < -M_PI) phil += m_2pi;
00160         philast = phil;
00161         if (rl > gvout) notout = 0;
00162         if (phil < phim) phim = phil;
00163         if (phil > phip) phip = phil;
00164       }
00165     }
00166   } else if ((rmin>gvin) && (rmax>gvout)) {
00167     if(5 == debug) std::cout<<__FILE__<<" case B (albedo) "<<std::endl;
00168     // case B (albedo)
00169     double len=0; double lstep=m_2pi*rhel/16.; if (lstep>10.)lstep=10.;
00170     double xl=temphel->Xh(len); double yl=temphel->Yh(len);
00171     double phil=atan2(yl,xl); double rl=sqrt(xl*xl+yl*yl); 
00172     phim=phil; phip=phil; double phis=phil; double deltap,dp1,dp2;
00173     int nstep=0; double philast=phil;
00174     while ((rl<gvout)&&(nstep<20)){
00175       len+=lstep; nstep++; xl=temphel->Xh(len); yl=temphel->Yh(len); 
00176       phil=atan2(yl,xl); rl=sqrt(xl*xl+yl*yl);
00177       deltap=phil-philast; if (deltap>M_PI)phil-=m_2pi;
00178       if (deltap<-M_PI)phil+=m_2pi; philast=phil;
00179       if (phil<phim)phim=phil; if (phil>phip)phip=phil;
00180     }
00181     dp1=fabs(phis-phim); dp2=fabs(phip-phis); deltap=dp1; 
00182     if(dp2>dp1)deltap=dp2; phim=phis-deltap; phip=phis+deltap;
00183   } else if ((rmin<gvin) && (rmax<gvout)) {
00184     // case C (looper)
00185     if(5 == debug) std::cout<<__FILE__<<" case C (looper) "<<std::endl;
00186     if (rc>rhel){
00187       double alp=asin(rhel/rc); phim=pc-alp; phip=pc+alp;
00188     }else{
00189       // going to have to step to find min, max
00190       double len=0; double lstep=m_2pi*rhel/16.; if (lstep>10.)lstep=10.;
00191       double xl=temphel->Xh(len); double yl=temphel->Yh(len);
00192       double phil=atan2(yl,xl); double rl=sqrt(xl*xl+yl*yl);
00193       int nstep=0; double philast=phil; int isin=0; int notout=1;
00194       while ((notout)&&(nstep<33)){
00195         len+=lstep; nstep++; xl=temphel->Xh(len); yl=temphel->Yh(len);
00196         phil=atan2(yl,xl); rl=sqrt(xl*xl+yl*yl);
00197         if ((rl<gvin)&&(!isin)){
00198           philast=phil;
00199         }else if((rl>gvin)&&(!isin)){
00200           isin=1; phim=philast; phip=philast;
00201         }
00202         if (isin){
00203           double deltap=phil-philast; if (deltap>M_PI)phil-=m_2pi;
00204           if (deltap<-M_PI)phil+=m_2pi; philast=phil; if (rl<gvin)notout=0;
00205           if (phil<phim)phim=phil; if (phil>phip)phip=phil;
00206         }
00207       }
00208     }//end rc<>rhel
00209   } else if ((rmin>gvin) && (rmax<gvout)) {
00210     // case D (contained)
00211     if (rc > rhel) {
00212       double alp = asin(rhel/rc);
00213       phim = pc - alp;
00214       phip = pc + alp;
00215     }
00216     // if rc<rhel and case D, it's an orbiter. Use max phim, phip, so do nothing
00217   }
00218   phim -= phiex;
00219   phip += phiex; 
00220 
00221   if(5 == debug) { 
00222     std::cout<<__FILE__<<" phim "<<phim <<" phip "<<phip<<  std::endl;
00223   }
00224 
00225   //  now try to add hits
00226   kkl = 0; 
00227   while(hhhh[kkl]) {
00228     MdcxHit* temphit = hhhh[kkl++];
00229     double factor=1.;
00230     if ((0 != temphit->type()) && (temphit->Layer()>7)) factor = MdcxParameters::addHitFactor;//yzhang FIXME 2009-11-03 
00231     if (5 == debug){
00232       std::cout<<__FILE__<<" "<<__LINE__<< " try to add hit:"<< std::endl;
00233       temphit->print(std::cout);
00234       std::cout<<__FILE__
00235         <<" used "<<temphit->GetUsedOnHel()
00236         <<" pull " << fabs(temphit->pull(*temphel))
00237         <<" <? nSig " << MdcxParameters::nSigAddHitTrk
00238         <<  std::endl;
00239     }
00240     //yzhang 2009-10-21 16:55:25  ///FIXME
00241     if( (!temphit->GetUsedOnHel()) 
00242         && (fabs(temphit->pull(*temphel))< factor * MdcxParameters::nSigAddHitTrk) ) {
00243       //if( (!temphit->GetUsedOnHel()) && (fabs(temphit->d())<1.2) ) //delete
00244       double pw = temphit->pw(); 
00245       if((phip-pw) > m_2pi)  pw += m_2pi;
00246       if((phim-pw) < -m_2pi) pw -= m_2pi;
00247       if ( (pw>phim) && (pw<phip) ) {
00248         temphit->SetConstErr(0);
00249         double pull = temphit->pull( *temphel );
00250         ncalc++;
00251         sumpull += fabs(pull);
00252         //cout << "MdcxAddHits: hit " << kkl-1 << " trk " << whichtrack << " pull " << pull << endl;
00253         if(g_addHitCut) g_addHitCut->fill(fabs(pull));
00254         int layer = temphit->Layer();
00255         if(g_addHitCut2d) g_addHitCut2d->fill(layer,fabs(pull));
00256         if(5 == debug) 
00257           std::cout<<__FILE__<<" "
00258             << " pull "<<pull 
00259             << " addcut "<< MdcxParameters::nSigAddHitTrk 
00260             << " len "<< temphel->Doca_Len()
00261             << " >? lmax "<< lmax
00262             << " >? maxTkLen "<< MdcxParameters::maxTrkLength
00263             <<  std::endl;
00264         //  yzhang 2009-10-21 22:55:26
00265         if(fabs(pull) < factor * MdcxParameters::nSigAddHitTrk) {
00266         //if(fabs(pull) < addcut) 
00267           double len = temphel->Doca_Len(); 
00268           //{ hhhh[kkl-1]->print(std::cout);
00269           //  cout << " trk " << whichtrack << " pull " << pull 
00270           //     << " d() " << temphit->d() << " len " << len << endl; }
00271           //int theflag = temphel->GetTurnFlag();  //zoujh ???
00272           if (((len>0.0)&&(len<lmax))||((len<0.0)&&(len>-MdcxParameters::maxTrkLength))) { 
00273             //templist.append(temphit); 
00274             listoasses.append(temphit);
00275             temphit->SetUsedOnHel(1);
00276             nadded++;
00277             if(5 == debug) std::cout<<"Added "<<  std::endl;
00278           }
00279         }
00280         temphit->SetConstErr(1);
00281       }//end phim,phip cuts
00282     } 
00283   }
00284   //listoasses = templist;
00285   return listoasses;
00286 }

Member Data Documentation

float MdcxAddHits::addcut [protected]
 

HepAList<MdcxHit> MdcxAddHits::hhhh [protected]
 

HepAList<MdcxHit> MdcxAddHits::hhhh [protected]
 

int MdcxAddHits::kkk [protected]
 

int MdcxAddHits::kkl [protected]
 

HepAList<MdcxHit> MdcxAddHits::listoasses [protected]
 

HepAList<MdcxHit> MdcxAddHits::listoasses [protected]
 

int MdcxAddHits::nadded [protected]
 

int MdcxAddHits::ncalc [protected]
 

int MdcxAddHits::nhits [protected]
 

int MdcxAddHits::ntracks [protected]
 

double MdcxAddHits::sumpull [protected]
 

int MdcxAddHits::thot [protected]
 

HepAList<MdcxHel> MdcxAddHits::tttt [protected]
 

HepAList<MdcxHel> MdcxAddHits::tttt [protected]
 

int MdcxAddHits::tuot [protected]
 

The documentation for this class was generated from the following files:
Generated on Wed Feb 2 16:30:04 2011 for BOSS6.5.5 by  doxygen 1.3.9.1