/[PAMELA software]/DarthVader/TrackerLevel2/src/TrkLevel2.cpp
ViewVC logotype

Diff of /DarthVader/TrackerLevel2/src/TrkLevel2.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.37 by pam-fi, Wed Jun 6 09:36:07 2007 UTC revision 1.42 by pam-fi, Tue Nov 27 15:28:57 2007 UTC
# Line 52  TrkTrack::TrkTrack(){ Line 52  TrkTrack::TrkTrack(){
52          ayv[ip]    = 0;          ayv[ip]    = 0;
53          dedx_x[ip] = 0;          dedx_x[ip] = 0;
54          dedx_y[ip] = 0;          dedx_y[ip] = 0;
55            multmaxx[ip] = 0;
56            multmaxy[ip] = 0;
57            seedx[ip]  = 0;  
58            seedy[ip]  = 0;
59            xpu[ip]    = 0;  
60            ypu[ip]    = 0;  
61    
62      };      };
63  //    clx = 0;  
64  //    cly = 0;  //     TrkParams::SetTrackingMode();
65  //    clx = new TRefArray(6,0); //forse causa memory leak???  //     TrkParams::SetPrecisionFactor();
66  //    cly = new TRefArray(6,0); //forse causa memory leak???  //     TrkParams::SetStepMin();
67  //    clx = TRefArray(6,0);      TrkParams::SetMiniDefault();
 //    cly = TRefArray(6,0);  
   
     TrkParams::SetTrackingMode();  
     TrkParams::SetPrecisionFactor();  
     TrkParams::SetStepMin();  
68      TrkParams::SetPFA();      TrkParams::SetPFA();
69    
70  };  };
# Line 96  TrkTrack::TrkTrack(const TrkTrack& t){ Line 98  TrkTrack::TrkTrack(const TrkTrack& t){
98          ayv[ip]    = t.ayv[ip];          ayv[ip]    = t.ayv[ip];
99          dedx_x[ip] = t.dedx_x[ip];          dedx_x[ip] = t.dedx_x[ip];
100          dedx_y[ip] = t.dedx_y[ip];          dedx_y[ip] = t.dedx_y[ip];
101            multmaxx[ip] = t.multmaxx[ip];
102            multmaxy[ip] = t.multmaxy[ip];
103            seedx[ip]    = t.seedx[ip];  
104            seedy[ip]    = t.seedy[ip];
105            xpu[ip]      = t.xpu[ip];  
106            ypu[ip]      = t.ypu[ip];  
107      };      };
108  //    clx = 0;  
109  //    cly = 0;  //     TrkParams::SetTrackingMode();
110  //    if(t.clx)clx = new TRefArray(*(t.clx));  //     TrkParams::SetPrecisionFactor();
111  //    if(t.cly)cly = new TRefArray(*(t.cly));  //     TrkParams::SetStepMin();  
112  //    clx = TRefArray(t.clx);      TrkParams::SetMiniDefault();
 //    cly = TRefArray(t.cly);  
   
     TrkParams::SetTrackingMode();  
     TrkParams::SetPrecisionFactor();  
     TrkParams::SetStepMin();    
113      TrkParams::SetPFA();      TrkParams::SetPFA();
114    
115  };  };
# Line 141  void TrkTrack::Copy(TrkTrack& t){ Line 144  void TrkTrack::Copy(TrkTrack& t){
144          t.ayv[ip]    = ayv[ip];          t.ayv[ip]    = ayv[ip];
145          t.dedx_x[ip] = dedx_x[ip];          t.dedx_x[ip] = dedx_x[ip];
146          t.dedx_y[ip] = dedx_y[ip];          t.dedx_y[ip] = dedx_y[ip];
147            t.multmaxx[ip] = multmaxx[ip];
148            t.multmaxy[ip] = multmaxy[ip];
149            t.seedx[ip]    = seedx[ip];  
150            t.seedy[ip]    = seedy[ip];
151            t.xpu[ip]      = xpu[ip];  
152            t.ypu[ip]      = ypu[ip];  
153                            
154      };      };
155    
 //    t.clx = TRefArray(clx);  
 //    t.cly = TRefArray(cly);  
156            
157  };  };
158  //--------------------------------------  //--------------------------------------
# Line 397  Float_t TrkTrack::GetLnLX(){ Line 404  Float_t TrkTrack::GetLnLX(){
404              lnl += (tailx[ip]+1.) * log( (tailx[ip]*pow(resx[ip],2.) + pow(xv[ip]-xm[ip],2.)) / (tailx[ip]*pow(resx[ip],2)) );              lnl += (tailx[ip]+1.) * log( (tailx[ip]*pow(resx[ip],2.) + pow(xv[ip]-xm[ip],2.)) / (tailx[ip]*pow(resx[ip],2)) );
405      if(GetNX()>3)lnl=lnl/(GetNX()-3);      if(GetNX()>3)lnl=lnl/(GetNX()-3);
406      else lnl=0;      else lnl=0;
407      if(lnl==0)cout << " Float_t  TrkTrack::GetLnLX() -- WARNING -- value not defined "<<lnl<<endl;      if(lnl==0){
408            cout << " Float_t  TrkTrack::GetLnLX() -- WARNING -- value not defined "<<lnl<<endl;
409            Dump();
410        }
411      return lnl;      return lnl;
412            
413  }  }
# Line 411  Float_t TrkTrack::GetLnLY(){ Line 421  Float_t TrkTrack::GetLnLY(){
421              lnl += (taily[ip]+1.) * log( (taily[ip]*pow(resy[ip],2.) + pow(yv[ip]-ym[ip],2.)) / (taily[ip]*pow(resy[ip],2)) );              lnl += (taily[ip]+1.) * log( (taily[ip]*pow(resy[ip],2.) + pow(yv[ip]-ym[ip],2.)) / (taily[ip]*pow(resy[ip],2)) );
422      if(GetNY()>2)lnl=lnl/(GetNY()-2);      if(GetNY()>2)lnl=lnl/(GetNY()-2);
423      else lnl=0;      else lnl=0;
424      if(lnl==0)cout << " Float_t  TrkTrack::GetLnLY() -- WARNING -- value not defined "<<lnl<<endl;      if(lnl==0){
425            cout << " Float_t  TrkTrack::GetLnLY() -- WARNING -- value not defined "<<lnl<<endl;
426            Dump();
427        }
428      return lnl;      return lnl;
429            
430  }  }
431    /**
432     * Returns the effective angle, relative to the sensor, on each plane.
433     * @param ip plane (0-5)
434     * @param iv view (0=x 1=y)
435     */
436    Float_t TrkTrack::GetEffectiveAngle(int ip, int iv){
437    
438        if(ip<0 || ip>5){
439            cout << "Float_t TrkTrack::GetEffectiveAngle(int "<<ip<<", int "<<iv<<") ==> wrong input"<<endl;
440            return 0.;
441        }
442    
443        float v[3]={xv[ip],yv[ip],zv[ip]};
444        //-----------------------------------------
445        // effective angle (relative to the sensor)
446        //-----------------------------------------
447        float axv_geo  = axv[ip];
448        float muhall_h = 297.61; //cm**2/Vs
449        float BY = TrkParams::GetBY(v);
450        float axv_eff = 0;
451        if(ip==5) axv_geo = -1*axv_geo;
452        if(ip==5) BY      = -1*BY;
453        axv_eff = 180.*atan( tan(axv_geo*acos(-1.)/180.) + muhall_h * BY * 0.0001)/acos(-1.);
454        //-----------------------------------------
455        // effective angle (relative to the sensor)
456        //-----------------------------------------
457        float ayv_geo = ayv[ip];
458        float muhall_e = 1258.18; //cm**2/Vs
459        float BX = TrkParams::GetBX(v);
460        float ayv_eff = 0;
461        ayv_eff = 180.*atan( tan(ayv_geo*acos(-1.)/180.) + muhall_e * BX * 0.0001)/acos(-1.);
462      
463        if     (iv==0)return axv_eff;
464        else if(iv==1)return ayv_eff;
465        else{
466            cout << "Float_t TrkTrack::GetEffectiveAngle(int "<<ip<<", int "<<iv<<") ==> wrong input"<<endl;
467            return 0.;
468        }
469      
470    };
471    
472  //--------------------------------------  //--------------------------------------
473  //  //
474  //  //
# Line 445  void TrkTrack::Dump(){ Line 499  void TrkTrack::Dump(){
499      cout << endl << "           "; for(int i=0; i<5; i++)cout << coval[4][i]<<" ";      cout << endl << "           "; for(int i=0; i<5; i++)cout << coval[4][i]<<" ";
500      cout << endl << "dedx_x   : "; for(int i=0; i<6; i++)cout << dedx_x[i] << " ";      cout << endl << "dedx_x   : "; for(int i=0; i<6; i++)cout << dedx_x[i] << " ";
501      cout << endl << "dedx_y   : "; for(int i=0; i<6; i++)cout << dedx_y[i] << " ";      cout << endl << "dedx_y   : "; for(int i=0; i<6; i++)cout << dedx_y[i] << " ";
502        cout << endl << "maxs x   : "; for(int i=0; i<6; i++)cout << GetClusterX_MaxStrip(i) << " ";
503        cout << endl << "maxs y   : "; for(int i=0; i<6; i++)cout << GetClusterY_MaxStrip(i) << " ";
504        cout << endl << "mult x   : "; for(int i=0; i<6; i++)cout << GetClusterX_Multiplicity(i) << " ";
505        cout << endl << "mult y   : "; for(int i=0; i<6; i++)cout << GetClusterY_Multiplicity(i) << " ";
506        cout << endl << "seed x   : "; for(int i=0; i<6; i++)cout << GetClusterX_Seed(i) << " ";
507        cout << endl << "seed y   : "; for(int i=0; i<6; i++)cout << GetClusterY_Seed(i) << " ";
508        cout << endl << "xpu      : "; for(int i=0; i<6; i++)cout << xpu[i] << " ";
509        cout << endl << "ypu      : "; for(int i=0; i<6; i++)cout << ypu[i] << " ";
510    
511      cout << endl;      cout << endl;
512  }  }
513  /**  /**
# Line 560  void TrkTrack::LoadField(TString path){ Line 623  void TrkTrack::LoadField(TString path){
623  //     path_.error   = 0;  //     path_.error   = 0;
624  //     readb_();  //     readb_();
625    
626      TrkParams::SetTrackingMode();  //     TrkParams::SetTrackingMode();
627      TrkParams::SetPrecisionFactor();  //     TrkParams::SetPrecisionFactor();
628      TrkParams::SetStepMin();  //     TrkParams::SetStepMin();
629        TrkParams::SetMiniDefault();
630    
631      TrkParams::Set(path,1);      TrkParams::Set(path,1);
632      TrkParams::Load(1);      TrkParams::Load(1);
# Line 663  void TrkTrack::SetFromMiniStruct(cMini2t Line 727  void TrkTrack::SetFromMiniStruct(cMini2t
727   */   */
728  Bool_t TrkTrack::EvaluateClusterPositions(){  Bool_t TrkTrack::EvaluateClusterPositions(){
729            
730  //     cout << "void TrkTrack::GetClusterPositions() "<<endl;  //     cout << "void TrkTrack::GetClusterositions() "<<endl;
731    
732      TrkParams::Load( );      TrkParams::Load( );
733      if( !TrkParams::IsLoaded() )return false;      if( !TrkParams::IsLoaded() )return false;
# Line 711  Bool_t TrkTrack::EvaluateClusterPosition Line 775  Bool_t TrkTrack::EvaluateClusterPosition
775   * @see EvaluateClusterPositions()   * @see EvaluateClusterPositions()
776   *   *
777   * The fitting procedure can be varied by changing the tracking mode,   * The fitting procedure can be varied by changing the tracking mode,
778   * the fit-precision factor and the minimum number of step.   * the fit-precision factor, the minimum number of step, etc.
779   * @see SetTrackingMode(int)   * @see SetTrackingMode(int)
780   * @see SetPrecisionFactor(double)   * @see SetPrecisionFactor(double)
781   * @see SetStepMin(int)   * @see SetStepMin(int)
782     * @see SetDeltaB(int,double)
783   */   */
784  void TrkTrack::Fit(double pfixed, int& fail, int iprint, int froml1){  void TrkTrack::Fit(double pfixed, int& fail, int iprint, int froml1){
785    
# Line 817  void TrkTrack::SetStepMin(int istepmin){ Line 882  void TrkTrack::SetStepMin(int istepmin){
882      track_.istepmin = istepmin;      track_.istepmin = istepmin;
883  }  }
884  /**  /**
885     * Set deltaB parameters (id=0,1). By default they are set to zero.
886     */
887    void TrkTrack::SetDeltaB(int id, double db){
888        if(id!=0 && id!=1)cout << "void TrkTrack::SetDeltaB(int id,double db) -- wrong input parameters: "<<id<<" "<<db<<endl;
889        TrkParams::SetDeltaB(id,db);
890    }
891    
892    /**
893   * Returns 1 if the track is inside the magnet cavity   * Returns 1 if the track is inside the magnet cavity
894   * Set the minimum number of steps for tracking precision   * Set the minimum number of steps for tracking precision
895   */   */
# Line 848  Int_t TrkTrack::GetClusterX_ID(int ip){ Line 921  Int_t TrkTrack::GetClusterX_ID(int ip){
921  Int_t TrkTrack::GetClusterY_ID(int ip){  Int_t TrkTrack::GetClusterY_ID(int ip){
922      return ((Int_t)fabs(ygood[ip]))%10000000-1;      return ((Int_t)fabs(ygood[ip]))%10000000-1;
923  };  };
924    
925  /**  /**
926   * Method to retrieve the ladder (0-4, increasing x) traversed by the track on this plane.   * Method to retrieve the ladder (0-4, increasing x) traversed by the track on this plane.
927   * If no ladder is traversed (dead area) the metod retuns -1.   * If no ladder is traversed (dead area) the metod retuns -1.
# Line 1031  void TrkLevel2::Dump(){ Line 1105  void TrkLevel2::Dump(){
1105                    
1106          //          //
1107      cout << endl << endl << "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-";      cout << endl << endl << "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-";
1108      cout << endl << "good     : "; for(int i=0; i<12; i++) cout << good[i]<<" ";      cout << endl << "good     : "; for(int i=0; i<12; i++) cout << hex <<" 0x"<< good[i]<<dec;
1109      cout << endl << "ntrk()   : " << this->ntrk() ;      cout << endl << "ntrk()   : " << ntrk() ;
1110      cout << endl << "nclsx()  : " << this->nclsx();      cout << endl << "nclsx()  : " << nclsx();
1111      cout << endl << "nclsy()  : " << this->nclsy();      cout << endl << "nclsy()  : " << nclsy();
1112      if(Track){      if(Track){
1113          TClonesArray &t  = *Track;          TClonesArray &t  = *Track;
1114          for(int i=0; i<ntrk(); i++)     ((TrkTrack *)t[i])->Dump();          for(int i=0; i<ntrk(); i++)     ((TrkTrack *)t[i])->Dump();
1115      }            }      
1116      if(SingletX){  //     if(SingletX){
1117          TClonesArray &sx = *SingletX;  //      TClonesArray &sx = *SingletX;
1118          for(int i=0; i<nclsx(); i++) ((TrkSinglet *)sx[i])->Dump();  //      for(int i=0; i<nclsx(); i++) ((TrkSinglet *)sx[i])->Dump();
1119      }  //     }
1120      if(SingletY){  //     if(SingletY){
1121          TClonesArray &sy = *SingletY;  //      TClonesArray &sy = *SingletY;
1122          for(int i=0; i<nclsy(); i++) ((TrkSinglet *)sy[i])->Dump();  //      for(int i=0; i<nclsy(); i++) ((TrkSinglet *)sy[i])->Dump();
1123      }  //     }
1124        cout << endl;
1125  }  }
1126  /**  /**
1127   * \brief Dump processing status   * \brief Dump processing status
# Line 1190  void TrkLevel2::SetFromLevel2Struct(cTrk Line 1265  void TrkLevel2::SetFromLevel2Struct(cTrk
1265              t_track->ayv[ip]    = l2->ayv_nt[i][ip];              t_track->ayv[ip]    = l2->ayv_nt[i][ip];
1266              t_track->dedx_x[ip] = l2->dedx_x[i][ip];              t_track->dedx_x[ip] = l2->dedx_x[i][ip];
1267              t_track->dedx_y[ip] = l2->dedx_y[i][ip];              t_track->dedx_y[ip] = l2->dedx_y[i][ip];
1268                t_track->multmaxx[ip] = l2->multmaxx[i][ip];
1269                t_track->multmaxy[ip] = l2->multmaxy[i][ip];
1270                t_track->seedx[ip]  = l2->seedx[i][ip];  
1271                t_track->seedy[ip]  = l2->seedy[i][ip];
1272                t_track->xpu[ip]    = l2->xpu[i][ip];  
1273                t_track->ypu[ip]    = l2->ypu[i][ip];  
1274              //-----------------------------------------------------              //-----------------------------------------------------
1275              //-----------------------------------------------------              //-----------------------------------------------------
1276              //-----------------------------------------------------              //-----------------------------------------------------
# Line 1404  TRefArray *TrkLevel2::GetTracks_NFitSort Line 1485  TRefArray *TrkLevel2::GetTracks_NFitSort
1485  TrkTrack *TrkLevel2::GetStoredTrack(int is){  TrkTrack *TrkLevel2::GetStoredTrack(int is){
1486    
1487      if(is >= this->ntrk()){      if(is >= this->ntrk()){
1488          cout << "** TrkLevel2 ** Track "<< is << "doen not exits! " << endl;          cout << "TrkTrack *TrkLevel2::GetStoredTrack(int) >> Track "<< is << "doen not exits! " << endl;
1489          cout << "                Stored tracks ntrk() = "<< this->ntrk() << endl;          cout << "Stored tracks ntrk() = "<< this->ntrk() << endl;
1490          return 0;          return 0;
1491      }      }
1492      if(!Track){      if(!Track){
# Line 1426  TrkTrack *TrkLevel2::GetStoredTrack(int Line 1507  TrkTrack *TrkLevel2::GetStoredTrack(int
1507  TrkSinglet *TrkLevel2::GetSingletX(int is){  TrkSinglet *TrkLevel2::GetSingletX(int is){
1508    
1509          if(is >= this->nclsx()){          if(is >= this->nclsx()){
1510                  cout << "** TrkLevel2 ** Singlet "<< is << "doen not exits! " << endl;                  cout << "TrkSinglet *TrkLevel2::GetSingletX(int) >> Singlet "<< is << "doen not exits! " << endl;
1511                  cout << "                Stored x-singlets nclsx() = "<< this->nclsx() << endl;                  cout << "Stored x-singlets nclsx() = "<< this->nclsx() << endl;
1512                  return 0;                  return 0;
1513          }          }
1514          if(!SingletX)return 0;          if(!SingletX)return 0;
# Line 1446  TrkSinglet *TrkLevel2::GetSingletX(int i Line 1527  TrkSinglet *TrkLevel2::GetSingletX(int i
1527  TrkSinglet *TrkLevel2::GetSingletY(int is){  TrkSinglet *TrkLevel2::GetSingletY(int is){
1528    
1529          if(is >= this->nclsy()){          if(is >= this->nclsy()){
1530                  cout << "** TrkLevel2 ** Singlet "<< is << "doen not exits! " << endl;                  cout << "TrkSinglet *TrkLevel2::GetSingletY(int) >> Singlet "<< is << "doen not exits! " << endl;
1531                  cout << "                Stored y-singlets nclsy() = "<< this->nclsx() << endl;                  cout << "Stored y-singlets nclsx() = "<< this->nclsx() << endl;
1532                  return 0;                  return 0;
1533          }          }
1534          if(!SingletY)return 0;          if(!SingletY)return 0;
# Line 1467  TrkSinglet *TrkLevel2::GetSingletY(int i Line 1548  TrkSinglet *TrkLevel2::GetSingletY(int i
1548  TrkTrack *TrkLevel2::GetTrack(int it){  TrkTrack *TrkLevel2::GetTrack(int it){
1549            
1550          if(it >= this->GetNTracks()){          if(it >= this->GetNTracks()){
1551                  cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl;                  cout << "TrkTrack *TrkLevel2::GetTrack(int) >> Track "<< it << "does not exits! " << endl;
1552                  cout << "                Physical tracks GetNTracks() = "<< this->ntrk() << endl;                  cout << "Physical tracks GetNTracks() = "<< this->ntrk() << endl;
1553                  return 0;                  return 0;
1554          }          }
1555                    
# Line 1505  Int_t TrkLevel2::GetNTracks(){ Line 1586  Int_t TrkLevel2::GetNTracks(){
1586  TrkTrack *TrkLevel2::GetTrackImage(int it){  TrkTrack *TrkLevel2::GetTrackImage(int it){
1587    
1588      if(it >= this->GetNTracks()){      if(it >= this->GetNTracks()){
1589          cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl;          cout << "TrkTrack *TrkLevel2::GetTrackImage(int) >> Track "<< it << "does not exits! " << endl;
1590          cout << "                Physical tracks GetNTracks() = "<< this->ntrk() << endl;          cout << "Physical tracks GetNTracks() = "<< this->ntrk() << endl;
1591          return 0;          return 0;
1592      }      }
1593                    
# Line 1515  TrkTrack *TrkLevel2::GetTrackImage(int i Line 1596  TrkTrack *TrkLevel2::GetTrackImage(int i
1596      TrkTrack *track = (TrkTrack*)sorted->At(it);      TrkTrack *track = (TrkTrack*)sorted->At(it);
1597                    
1598      if(!track->HasImage()){      if(!track->HasImage()){
1599          cout << "** TrkLevel2 ** Track "<< it << "does not have image! " << endl;          cout << "TrkTrack *TrkLevel2::GetTrackImage(int) >> Track "<< it << "does not have image! " << endl;
1600          return 0;          return 0;
1601      }      }
1602      if(!Track)return 0;      if(!Track)return 0;
# Line 1542  void TrkLevel2::LoadField(TString path){ Line 1623  void TrkLevel2::LoadField(TString path){
1623  //     path_.error   = 0;  //     path_.error   = 0;
1624  //     readb_();  //     readb_();
1625    
1626      TrkParams::SetTrackingMode();  //     TrkParams::SetTrackingMode();
1627      TrkParams::SetPrecisionFactor();  //     TrkParams::SetPrecisionFactor();
1628      TrkParams::SetStepMin();  //     TrkParams::SetStepMin();
1629        TrkParams::SetMiniDefault();
1630    
1631      TrkParams::Set(path,1);      TrkParams::Set(path,1);
1632      TrkParams::Load(1);      TrkParams::Load(1);
Legend:
Removed from v.1.37  
changed lines
  Added in v.1.42
  ViewVC Help
Powered by ViewVC 1.1.23