/[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.54 by pam-fi, Wed Mar 11 14:19:10 2009 UTC revision 1.58 by mocchiut, Wed Aug 27 14:44:37 2014 UTC
# Line 391  Int_t TrkTrack::GetLeverArmXY(){ Line 391  Int_t TrkTrack::GetLeverArmXY(){
391      return (last_plane-first_plane+1);      return (last_plane-first_plane+1);
392  }  }
393  /**  /**
394     * Returns the number of hit planes
395     */
396    Int_t TrkTrack::GetNhit()  {
397      int np=0;
398      for(Int_t ip=0; ip<6; ip++) np += (XGood(ip)||YGood(ip)) ;
399      return np;
400    };
401    /**
402   * Returns the reduced chi-square of track x-projection   * Returns the reduced chi-square of track x-projection
403   */   */
404  Float_t  TrkTrack::GetChi2X(){  Float_t  TrkTrack::GetChi2X(){
# Line 398  Float_t  TrkTrack::GetChi2X(){ Line 406  Float_t  TrkTrack::GetChi2X(){
406      for(int ip=0; ip<6; ip++)if(XGood(ip))chiq+= pow((xv[ip]-xm[ip])/resx[ip],2.);      for(int ip=0; ip<6; ip++)if(XGood(ip))chiq+= pow((xv[ip]-xm[ip])/resx[ip],2.);
407      if(GetNX()>3)chiq=chiq/(GetNX()-3);      if(GetNX()>3)chiq=chiq/(GetNX()-3);
408      else chiq=0;      else chiq=0;
409      if(chiq==0)cout << " Float_t  TrkTrack::GetChi2X() -- WARNING -- value not defined "<<chiq<<endl;      //    if(chiq==0)cout << " Float_t  TrkTrack::GetChi2X() -- WARNING -- value not defined "<<chiq<<endl;
410      return chiq;      return chiq;
411  }  }
412  /**  /**
# Line 409  Float_t  TrkTrack::GetChi2Y(){ Line 417  Float_t  TrkTrack::GetChi2Y(){
417      for(int ip=0; ip<6; ip++)if(YGood(ip))chiq+= pow((yv[ip]-ym[ip])/resy[ip],2.);      for(int ip=0; ip<6; ip++)if(YGood(ip))chiq+= pow((yv[ip]-ym[ip])/resy[ip],2.);
418      if(GetNY()>2)chiq=chiq/(GetNY()-2);      if(GetNY()>2)chiq=chiq/(GetNY()-2);
419      else chiq=0;      else chiq=0;
420      if(chiq==0)cout << " Float_t  TrkTrack::GetChi2Y() -- WARNING -- value not defined "<<chiq<<endl;      //    if(chiq==0)cout << " Float_t  TrkTrack::GetChi2Y() -- WARNING -- value not defined "<<chiq<<endl;
421      return chiq;      return chiq;
422  }  }
423  /**  /**
# Line 664  void TrkTrack::FillMiniStruct(cMini2trac Line 672  void TrkTrack::FillMiniStruct(cMini2trac
672    
673  //      cout << i<<" - "<<xgood[i]<<" "<<XGood(i)<<endl;  //      cout << i<<" - "<<xgood[i]<<" "<<XGood(i)<<endl;
674  //      cout << i<<" - "<<ygood[i]<<" "<<YGood(i)<<endl;  //      cout << i<<" - "<<ygood[i]<<" "<<YGood(i)<<endl;
675          track.xgood[i]=XGood(i);        track.xgood[i] = (double) XGood(i);
676          track.ygood[i]=YGood(i);          track.ygood[i] = (double) YGood(i);
677                    
678          track.xm[i]=xm[i];          track.xm[i]= (double) xm[i];
679          track.ym[i]=ym[i];          track.ym[i]= (double) ym[i];
680          track.zm[i]=zm[i];          track.zm[i]= (double) zm[i];
681                    
682  //      --- temporaneo ----------------------------  //      --- temporaneo ----------------------------
683  //      float segment = 100.;  //      float segment = 100.;
# Line 687  void TrkTrack::FillMiniStruct(cMini2trac Line 695  void TrkTrack::FillMiniStruct(cMini2trac
695  //      --- temporaneo ----------------------------  //      --- temporaneo ----------------------------
696    
697          if( XGood(i) || YGood(i) ){          if( XGood(i) || YGood(i) ){
698              double segment = 2.;//cm              //NB!! the length of the sensor is not exactely taken into account    
699                double segment = 7.;// 2.;//cm //Elena 10th
700              // NB: i parametri di allineamento hanno una notazione particolare!!!              // NB: i parametri di allineamento hanno una notazione particolare!!!
701              // sensor = 0 (hybrid side), 1              // sensor = 0 (hybrid side), 1
702              // ladder = 0-2 (increasing x)              // ladder = 0-2 (increasing x)
# Line 697  void TrkTrack::FillMiniStruct(cMini2trac Line 706  void TrkTrack::FillMiniStruct(cMini2trac
706              int il = (int)GetLadder(i);              int il = (int)GetLadder(i);
707                            
708              double omega   = 0.;              double omega   = 0.;
709              double beta    = 0.;              //      double beta    = 0.;// EM GCC 4.7
710              double gamma   = 0.;              //      double gamma   = 0.;
711              if(              if(
712                  (is < 0 || is > 1 || ip < 0 || ip > 5 || il < 0 || il > 2) &&                  (is < 0 || is > 1 || ip < 0 || ip > 5 || il < 0 || il > 2) &&
713                  true){                  true){
# Line 708  void TrkTrack::FillMiniStruct(cMini2trac Line 717  void TrkTrack::FillMiniStruct(cMini2trac
717                  cout << " is ip il = "<<is<<" "<<ip<<" "<<il<<endl;                  cout << " is ip il = "<<is<<" "<<ip<<" "<<il<<endl;
718              }else{              }else{
719                  omega   = alignparameters_.omega[is][il][ip];                  omega   = alignparameters_.omega[is][il][ip];
720                  beta    = alignparameters_.beta[is][il][ip];                  //              beta    = alignparameters_.beta[is][il][ip];// EM GCC 4.7 unused
721                  gamma   = alignparameters_.gamma[is][il][ip];                  //              gamma   = alignparameters_.gamma[is][il][ip];// EM GCC 4.7 unused
722              }              }
723                            
724              if(       XGood(i) && !YGood(i) ){              if(       XGood(i) && !YGood(i) ){
725                  track.xm_a[i] = xm[i] - omega * segment;                  track.xm_a[i] =  (double) xm[i] - omega * segment;
726                  track.ym_a[i] = ym[i] + segment;                  track.ym_a[i] =  (double) ym[i] + segment;
727  //          track.zm_a[i] = zm[i] + beta * segment;//not used yet  //          track.zm_a[i] = zm[i] + beta * segment;//not used yet
728                  track.xm_b[i] = xm[i] + omega * segment;                  track.xm_b[i] =  (double) xm[i] + omega * segment;
729                  track.ym_b[i] = ym[i] - segment;                  track.ym_b[i] =  (double) ym[i] - segment;
730  //          track.zm_b[i] = zm[i] - beta * segment;//not used yet  //          track.zm_b[i] = zm[i] - beta * segment;//not used yet
731              }else if( !XGood(i) && YGood(i) ){              }else if( !XGood(i) && YGood(i) ){
732                  track.xm_a[i] = xm[i] + segment;                  track.xm_a[i] =  (double) xm[i] + segment;
733                  track.ym_a[i] = ym[i] + omega * segment;                  track.ym_a[i] =  (double) ym[i] + omega * segment;
734  //          track.zm_a[i] = zm[i] - gamma * segment;//not used yet  //          track.zm_a[i] = zm[i] - gamma * segment;//not used yet
735                  track.xm_b[i] = xm[i] - segment;                  track.xm_b[i] =  (double) xm[i] - segment;
736                  track.ym_b[i] = ym[i] - omega * segment;                  track.ym_b[i] =  (double) ym[i] - omega * segment;
737  //          track.zm_b[i] = zm[i] + gamma * segment;//not used yet  //          track.zm_b[i] = zm[i] + gamma * segment;//not used yet
738              }              }
739          }          }
740                    
741          track.resx[i]=resx[i];          track.resx[i]= (double) resx[i];
742          track.resy[i]=resy[i];          track.resy[i]= (double) resy[i];
743          track.tailx[i]=tailx[i];          track.tailx[i]= (double) tailx[i];
744          track.taily[i]=taily[i];          track.taily[i]= (double) taily[i];
745      }      }
746    
747      for(int i=0; i<5; i++) track.al[i]=al[i];      for(int i=0; i<5; i++) track.al[i]= (double) al[i];
748      track.zini = 23.5;      track.zini = 23.5;
749  // ZINI = 23.5 !!! it should be the same parameter in all codes  // ZINI = 23.5 !!! it should be the same parameter in all codes
750            
# Line 746  void TrkTrack::FillMiniStruct(cMini2trac Line 755  void TrkTrack::FillMiniStruct(cMini2trac
755  void TrkTrack::SetFromMiniStruct(cMini2track *track){  void TrkTrack::SetFromMiniStruct(cMini2track *track){
756    
757      for(int i=0; i<5; i++) {      for(int i=0; i<5; i++) {
758          al[i]=track->al[i];        al[i]= (float) (track->al[i]);
759          for(int j=0; j<5; j++) coval[i][j]=track->cov[i][j];          for(int j=0; j<5; j++) coval[i][j]= (float) (track->cov[i][j]);
760      }      }
761      chi2  = track->chi2;      chi2  =  (float) (track->chi2);
762      nstep = track->nstep;      nstep =  (float) (track->nstep);
763      for(int i=0; i<6; i++){      for(int i=0; i<6; i++){
764          xv[i]  = track->xv[i];          xv[i]  =  (float) (track->xv[i]);
765          yv[i]  = track->yv[i];          yv[i]  =  (float) (track->yv[i]);
766          zv[i]  = track->zv[i];          zv[i]  =  (float) (track->zv[i]);
767          xm[i]  = track->xm[i];          xm[i]  = (float)  (track->xm[i]);
768          ym[i]  = track->ym[i];          ym[i]  =  (float) (track->ym[i]);
769          zm[i]  = track->zm[i];          zm[i]  =  (float) (track->zm[i]);
770          axv[i] = track->axv[i];          axv[i] =  (float) (track->axv[i]);
771          ayv[i] = track->ayv[i];          ayv[i] =  (float) (track->ayv[i]);      
772            resx[i] = (float)  (track->resx[i]); //Elena 10th
773            resy[i] =  (float) (track->resy[i]);
774      }      }
775            
776  }  }
# Line 810  Bool_t TrkTrack::EvaluateClusterPosition Line 821  Bool_t TrkTrack::EvaluateClusterPosition
821    
822      for(int ip=0; ip<6; ip++){      for(int ip=0; ip<6; ip++){
823  //      cout << ip<<" ** "<<xm[ip]<<" / "<<ym[ip]<<endl;;  //      cout << ip<<" ** "<<xm[ip]<<" / "<<ym[ip]<<endl;;
824          int icx = GetClusterX_ID(ip)+1;          int icx = GetClusterX_ID(ip)+1;//0=no-cluster,1-N
825          int icy = GetClusterY_ID(ip)+1;          int icy = GetClusterY_ID(ip)+1;//0=no-cluster,1-N
826          int sensor = GetSensor(ip)+1;//<< convenzione "Paolo"          int sensor = GetSensor(ip)+1;//<< convenzione "Paolo"
827          if(ip==5 && sensor!=0)sensor=3-sensor;//<< convenzione "Elena"          if(ip==5 && sensor!=0)sensor=3-sensor;//<< convenzione "Elena"
828          int ladder = GetLadder(ip)+1;          int ladder = GetLadder(ip)+1;
# Line 825  Bool_t TrkTrack::EvaluateClusterPosition Line 836  Bool_t TrkTrack::EvaluateClusterPosition
836          float bfy = 10*TrkParams::GetBY(v);//Tesla          float bfy = 10*TrkParams::GetBY(v);//Tesla
837          int ipp=ip+1;          int ipp=ip+1;
838          xyzpam_(&ipp,&icx,&icy,&ladder,&sensor,&ax,&ay,&bfx,&bfy);          xyzpam_(&ipp,&icx,&icy,&ladder,&sensor,&ax,&ay,&bfx,&bfy);
839          if(icx<0 || icy<0)return false;          //      if(icx<0 || icy<0)return false;
840      }      }
841      return true;      return true;
842  }  }
# Line 875  void TrkTrack::Fit(double pfixed, int& f Line 886  void TrkTrack::Fit(double pfixed, int& f
886      extern cMini2track track_;      extern cMini2track track_;
887      fail = 0;      fail = 0;
888    
889      FillMiniStruct(track_);      //    FillMiniStruct(track_);
890                    
891      if(froml1!=0){      if(froml1!=0){
892          if( !EvaluateClusterPositions() ){          if( !EvaluateClusterPositions() ){
# Line 914  void TrkTrack::Fit(double pfixed, int& f Line 925  void TrkTrack::Fit(double pfixed, int& f
925          if(iprint)cout << "ERROR: ifail= " << ifail << endl;          if(iprint)cout << "ERROR: ifail= " << ifail << endl;
926          fail = 1;          fail = 1;
927      }      }
928        if(chi2!=chi2){
929            if(iprint)cout << "ERROR: chi2= " << chi2 << endl;      
930            FitReset();
931            fail = 1;      
932        }
933      //  ------------------------------------------      //  ------------------------------------------
934            
935      SetFromMiniStruct(&track_);      SetFromMiniStruct(&track_);
# Line 1152  Int_t TrkTrack::GetSensor(int ip){ Line 1168  Int_t TrkTrack::GetSensor(int ip){
1168  void TrkTrack::SetXGood(int ip, int clid, int il, int is, bool bad){  void TrkTrack::SetXGood(int ip, int clid, int il, int is, bool bad){
1169  //    int il=0;       //ladder (temporary)  //    int il=0;       //ladder (temporary)
1170  //    bool bad=false; //ladder (temporary)  //    bool bad=false; //ladder (temporary)
1171      if(ip<0||ip>5||clid<0||il<-1||il>2||is<-1||is>1)      if(ip<0||ip>5||clid<1||il<-1||il>2||is<-1||is>1)
1172          cout << " void TrkTrack::SetXGood(int,int,int,int,bool) --> MA SEI DI COCCIO?!?!"<<endl;          cout << " void TrkTrack::SetXGood(int,int,int,int,bool) --> MA SEI DI COCCIO?!?!"<<endl;
1173      xgood[ip]=(il+1)*100000000+(is+1)*10000000+clid;      xgood[ip]=(il+1)*100000000+(is+1)*10000000+clid;
1174      if(bad)xgood[ip]=-xgood[ip];      if(bad)xgood[ip]=-xgood[ip];
# Line 1169  void TrkTrack::SetXGood(int ip, int clid Line 1185  void TrkTrack::SetXGood(int ip, int clid
1185  void TrkTrack::SetYGood(int ip, int clid, int il, int is, bool bad){  void TrkTrack::SetYGood(int ip, int clid, int il, int is, bool bad){
1186  //    int il=0;       //ladder (temporary)  //    int il=0;       //ladder (temporary)
1187  //    bool bad=false; //ladder (temporary)  //    bool bad=false; //ladder (temporary)
1188      if(ip<0||ip>5||clid<0||il<-1||il>2||is<-1||is>1)      if(ip<0||ip>5||clid<1||il<-1||il>2||is<-1||is>1)
1189          cout << " void TrkTrack::SetYGood(int,int,int,int,bool) --> MA SEI DI COCCIO?!?!"<<endl;          cout << " void TrkTrack::SetYGood(int,int,int,int,bool) --> MA SEI DI COCCIO?!?!"<<endl;
1190      ygood[ip]=(il+1)*100000000+(is+1)*10000000+clid;      ygood[ip]=(il+1)*100000000+(is+1)*10000000+clid;
1191      if(bad)ygood[ip]=-ygood[ip];      if(bad)ygood[ip]=-ygood[ip];
# Line 2178  Float_t TrkLevel2::GetZTrk(Int_t plane_i Line 2194  Float_t TrkLevel2::GetZTrk(Int_t plane_i
2194   * (By default is created with z-coordinates inside the tracking volume)   * (By default is created with z-coordinates inside the tracking volume)
2195    */    */
2196  Trajectory::Trajectory(){  Trajectory::Trajectory(){
2197      npoint = 10;      npoint = 6;
2198      x = new float[npoint];      x = new float[npoint];
2199      y = new float[npoint];      y = new float[npoint];
2200      z = new float[npoint];      z = new float[npoint];
# Line 2245  Trajectory::Trajectory(int n, float* zin Line 2261  Trajectory::Trajectory(int n, float* zin
2261      thy = new float[npoint];      thy = new float[npoint];
2262      tl = new float[npoint];      tl = new float[npoint];
2263      int i=0;      int i=0;
2264      do{      do{      
2265          x[i] = 0;          x[i] = 0.;
2266          y[i] = 0;          y[i] = 0.;
2267          z[i] = zin[i];          z[i] = zin[i];
2268          thx[i] = 0;          thx[i] = 0.;
2269          thy[i] = 0;          thy[i] = 0.;
2270          tl[i] = 0;          tl[i] = 0.;
2271          i++;                      i++;            
2272      }while(zin[i-1] > zin[i] && i < npoint);      }while(zin[i-1] > zin[i] && i < npoint);
2273      npoint=i;      npoint=i;
2274      if(npoint != n)cout << "NB! Trajectory created with "<<npoint<<" points"<<endl;      if(npoint != n)cout << "NB! Trajectory created with "<<npoint<<" points instean of "<<n<<endl;
2275        //    Dump();
2276  }  }
2277  void Trajectory::Delete(){  void Trajectory::Delete(){
2278            
Legend:
Removed from v.1.54  
changed lines
  Added in v.1.58
  ViewVC Help
Powered by ViewVC 1.1.23