/[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.21 by pam-fi, Thu Jan 11 10:20:58 2007 UTC revision 1.31 by pam-fi, Wed Mar 28 09:22:28 2007 UTC
# Line 13  extern "C" {     Line 13  extern "C" {    
13      void dotrack_(int*, double*, double*, double*, double*, int*);      void dotrack_(int*, double*, double*, double*, double*, int*);
14      void dotrack2_(int*, double*, double*, double*, double*,double*, double*, double*,int*);      void dotrack2_(int*, double*, double*, double*, double*,double*, double*, double*,int*);
15  //    int  readb_(const char*);  //    int  readb_(const char*);
16      int  readb_();  //    int  readb_();
17      void mini2_(int*,int*,int*);       void mini2_(int*,int*,int*);
18      void guess_();       void guess_();
19         void gufld_(float*, float*);
20  }  }
21    
22  //--------------------------------------  //--------------------------------------
23  //  //
24  //  //
# Line 49  TrkTrack::TrkTrack(){ Line 51  TrkTrack::TrkTrack(){
51      };      };
52      clx = 0;      clx = 0;
53      cly = 0;      cly = 0;
54    //    clx = new TRefArray(6,0); //forse causa memory leak???
55    //    cly = new TRefArray(6,0); //forse causa memory leak???
56  };  };
57  //--------------------------------------  //--------------------------------------
58  //  //
# Line 79  TrkTrack::TrkTrack(const TrkTrack& t){ Line 83  TrkTrack::TrkTrack(const TrkTrack& t){
83          dedx_x[ip] = t.dedx_x[ip];          dedx_x[ip] = t.dedx_x[ip];
84          dedx_y[ip] = t.dedx_y[ip];          dedx_y[ip] = t.dedx_y[ip];
85      };      };
86      clx = new TRefArray(*(t.clx));      clx = 0;
87      cly = new TRefArray(*(t.cly));      cly = 0;
88        if(t.clx)clx = new TRefArray(*(t.clx));
89        if(t.cly)cly = new TRefArray(*(t.cly));
90                    
91  };  };
92  //--------------------------------------  //--------------------------------------
# Line 139  int TrkTrack::DoTrack(Trajectory* t){ Line 145  int TrkTrack::DoTrack(Trajectory* t){
145      for (int i=0; i<5; i++)         dal[i]  = (double)al[i];      for (int i=0; i<5; i++)         dal[i]  = (double)al[i];
146      for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i];      for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i];
147    
148        TrkParams::Load(1);
149        if( !TrkParams::IsLoaded(1) ){
150            cout << "int TrkTrack::DoTrack(Trajectory* t) --- ERROR --- m.field not loaded"<<endl;
151            return 0;
152        }
153      dotrack_(&(t->npoint),dzin,dxout,dyout,dal,&ifail);      dotrack_(&(t->npoint),dzin,dxout,dyout,dal,&ifail);
154            
155      for (int i=0; i<t->npoint; i++){      for (int i=0; i<t->npoint; i++){
# Line 178  int TrkTrack::DoTrack2(Trajectory* t){ Line 189  int TrkTrack::DoTrack2(Trajectory* t){
189      for (int i=0; i<5; i++)         dal[i]  = (double)al[i];      for (int i=0; i<5; i++)         dal[i]  = (double)al[i];
190      for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i];      for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i];
191    
192        TrkParams::Load(1);
193        if( !TrkParams::IsLoaded(1) ){
194            cout << "int TrkTrack::DoTrack2(Trajectory* t) --- ERROR --- m.field not loaded"<<endl;
195            return 0;
196        }
197      dotrack2_(&(t->npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail);      dotrack2_(&(t->npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail);
198            
199      for (int i=0; i<t->npoint; i++){      for (int i=0; i<t->npoint; i++){
# Line 219  Float_t TrkTrack::GetDeflection(){ Line 235  Float_t TrkTrack::GetDeflection(){
235  //  //
236  Float_t TrkTrack::GetDEDX(){  Float_t TrkTrack::GetDEDX(){
237          Float_t dedx=0;          Float_t dedx=0;
238          for(Int_t i=0; i<6; i++)dedx+=dedx_x[i]*xgood[i]+dedx_y[i]*ygood[i];  //      for(Int_t i=0; i<6; i++)dedx+=dedx_x[i]*xgood[i]+dedx_y[i]*ygood[i];
239            for(Int_t i=0; i<6; i++)dedx+=dedx_x[i]*XGood(i)+dedx_y[i]*YGood(i);
240          dedx = dedx/(this->GetNX()+this->GetNY());          dedx = dedx/(this->GetNX()+this->GetNY());
241          return dedx;          return dedx;
242  };  };
# Line 235  void TrkTrack::Dump(){ Line 252  void TrkTrack::Dump(){
252      cout << endl << "al       : "; for(int i=0; i<5; i++)cout << al[i] << " ";      cout << endl << "al       : "; for(int i=0; i<5; i++)cout << al[i] << " ";
253      cout << endl << "chi^2    : "<< chi2;      cout << endl << "chi^2    : "<< chi2;
254      cout << endl << "n.step   : "<< nstep;      cout << endl << "n.step   : "<< nstep;
255      cout << endl << "xgood    : "; for(int i=0; i<6; i++)cout << xgood[i] ;      cout << endl << "xgood    : "; for(int i=0; i<6; i++)cout << XGood(i) ;
256      cout << endl << "ygood    : "; for(int i=0; i<6; i++)cout << ygood[i] ;      cout << endl << "ygood    : "; for(int i=0; i<6; i++)cout << YGood(i) ;
257      cout << endl << "xm       : "; for(int i=0; i<6; i++)cout << xm[i] << " ";      cout << endl << "xm       : "; for(int i=0; i<6; i++)cout << xm[i] << " ";
258      cout << endl << "ym       : "; for(int i=0; i<6; i++)cout << ym[i] << " ";      cout << endl << "ym       : "; for(int i=0; i<6; i++)cout << ym[i] << " ";
259      cout << endl << "zm       : "; for(int i=0; i<6; i++)cout << zm[i] << " ";      cout << endl << "zm       : "; for(int i=0; i<6; i++)cout << zm[i] << " ";
# Line 273  void TrkTrack::SetResolution(double *rx, Line 290  void TrkTrack::SetResolution(double *rx,
290   * Set the TrkTrack good measurement   * Set the TrkTrack good measurement
291   */   */
292  void TrkTrack::SetGood(int *xg, int *yg){  void TrkTrack::SetGood(int *xg, int *yg){
293        // NB! si perdera` l'informazione sul numero del cluster
294      for(int i=0; i<6; i++) xgood[i]=*xg++;      for(int i=0; i<6; i++) xgood[i]=*xg++;
295      for(int i=0; i<6; i++) ygood[i]=*yg++;      for(int i=0; i<6; i++) ygood[i]=*yg++;
296  }  }
# Line 282  void TrkTrack::SetGood(int *xg, int *yg) Line 300  void TrkTrack::SetGood(int *xg, int *yg)
300   */   */
301  void TrkTrack::LoadField(TString path){  void TrkTrack::LoadField(TString path){
302            
303      strcpy(path_.path,path.Data());  //     strcpy(path_.path,path.Data());
304      path_.pathlen = path.Length();  //     path_.pathlen = path.Length();
305      path_.error   = 0;  //     path_.error   = 0;
306      readb_();  //     readb_();
307    
308        TrkParams::Set(path,1);
309        TrkParams::Load(1);
310    
311  };  };
312    
313    
314  /**  /**
315   * Method to fill minimization-routine common   * Method to fill minimization-routine common
316   */   */
# Line 296  void TrkTrack::FillMiniStruct(cMini2trac Line 318  void TrkTrack::FillMiniStruct(cMini2trac
318    
319      for(int i=0; i<6; i++){      for(int i=0; i<6; i++){
320    
321          track.xgood[i]=xgood[i];  //      track.xgood[i]=xgood[i];
322          track.ygood[i]=ygood[i];  //      track.ygood[i]=ygood[i];
323            track.xgood[i]=XGood(i);
324            track.ygood[i]=YGood(i);
325                    
326          track.xm[i]=xm[i];          track.xm[i]=xm[i];
327          track.ym[i]=ym[i];          track.ym[i]=ym[i];
# Line 310  void TrkTrack::FillMiniStruct(cMini2trac Line 334  void TrkTrack::FillMiniStruct(cMini2trac
334          track.xm_b[i]=xm[i];          track.xm_b[i]=xm[i];
335          track.ym_a[i]=ym[i];          track.ym_a[i]=ym[i];
336          track.ym_b[i]=ym[i];          track.ym_b[i]=ym[i];
337          if(       xgood[i] && !ygood[i] ){          if(       XGood(i) && !YGood(i) ){
338              track.ym_a[i] = track.ym_a[i]+segment;              track.ym_a[i] = track.ym_a[i]+segment;
339              track.ym_b[i] = track.ym_b[i]-segment;              track.ym_b[i] = track.ym_b[i]-segment;
340          }else if( !xgood[i] && ygood[i]){          }else if( !XGood(i) && YGood(i)){
341              track.xm_a[i] = track.xm_a[i]+segment;              track.xm_a[i] = track.xm_a[i]+segment;
342              track.xm_b[i] = track.xm_b[i]-segment;              track.xm_b[i] = track.xm_b[i]-segment;
343          }          }
# Line 349  void TrkTrack::SetFromMiniStruct(cMini2t Line 373  void TrkTrack::SetFromMiniStruct(cMini2t
373          axv[i] = track->axv[i];          axv[i] = track->axv[i];
374          ayv[i] = track->ayv[i];          ayv[i] = track->ayv[i];
375      }      }
   
376            
377  }  }
378  /**  /**
# Line 381  void TrkTrack::Fit(double pfixed, int& f Line 404  void TrkTrack::Fit(double pfixed, int& f
404    
405      //  ------------------------------------------      //  ------------------------------------------
406      //  call mini routine      //  call mini routine
407        TrkParams::Load(1);
408        if( !TrkParams::IsLoaded(1) ){
409            cout << "void TrkTrack::Fit(double pfixed, int& fail, int iprint) --- ERROR --- m.field not loaded"<<endl;
410            return;
411        }
412      int istep=0;      int istep=0;
413      int ifail=0;      int ifail=0;
414      mini2_(&istep,&ifail, &iprint);      mini2_(&istep,&ifail, &iprint);
# Line 427  void TrkTrack::FitReset(){ Line 455  void TrkTrack::FitReset(){
455          for(int j=0; j<5; j++) coval[i][j]=0.;          for(int j=0; j<5; j++) coval[i][j]=0.;
456      }      }
457  }  }
458    /*
459     * Set the tracking mode
460     */
461    void TrkTrack::SetTrackingMode(int trackmode){
462        extern cMini2track track_;
463        track_.trackmode = trackmode;
464    }
465    /*
466     * Set the factor scale for tracking precision
467     */
468    void TrkTrack::SetPrecisionFactor(double fact){
469        extern cMini2track track_;
470        track_.fact = fact;
471    }
472    /*
473     * Set the factor scale for tracking precision
474     */
475    void TrkTrack::SetStepMin(int istepmin){
476        extern cMini2track track_;
477        track_.istepmin = istepmin;
478    }
479    
480    
481    /*
482     * Method to retrieve the X-view clusters associated to the track.
483     * @param ip Tracker plane (0-5)
484     */
485    TrkCluster *TrkTrack::GetClusterX(int ip){
486        cout << " TrkCluster *TrkTrack::GetClusterX(int ip) -- momentaneamente fuori servizio --"<< endl;
487        if(!clx)return NULL;
488        TrkCluster *pt = (TrkCluster*)(clx->At(ip));
489        return pt;
490    };
491    /*
492     * Method to retrieve the Y-view clusters associated to the track.
493     * @param ip Tracker plane (0-5)
494     */
495    TrkCluster *TrkTrack::GetClusterY(int ip){
496        cout << " TrkCluster *TrkTrack::GetClusterY(int ip) -- momentaneamente fuori servizio --"<< endl;
497        if(!cly)return NULL;
498        TrkCluster *pt = (TrkCluster*)(cly->At(ip));
499        return pt;
500    };
501    
502    
503  //--------------------------------------  //--------------------------------------
504  //  //
# Line 531  void TrkSinglet::Clear(){ Line 603  void TrkSinglet::Clear(){
603  //  //
604  //--------------------------------------  //--------------------------------------
605  TrkLevel2::TrkLevel2(){  TrkLevel2::TrkLevel2(){
606  //    cout <<"TrkLevel2::TrkLevel2()"<<endl;    //    cout <<"TrkLevel2::TrkLevel2()"<<endl;
607      for(Int_t i=0; i<12 ; i++){      for(Int_t i=0; i<12 ; i++){
608                  good[i] = -1;                  good[i] = -1;
609          };          };
610  //    Track    = new TClonesArray("TrkTrack");  // okkio!! memory-leak
611  //    SingletX = new TClonesArray("TrkSinglet");  //     Track    = new TClonesArray("TrkTrack");
612  //    SingletY = new TClonesArray("TrkSinglet");  //     SingletX = new TClonesArray("TrkSinglet");
613    //     SingletY = new TClonesArray("TrkSinglet");
614      Track    = 0;      Track    = 0;
615      SingletX = 0;      SingletX = 0;
616      SingletY = 0;      SingletY = 0;
# Line 547  TrkLevel2::TrkLevel2(){ Line 620  TrkLevel2::TrkLevel2(){
620  //  //
621  //  //
622  //--------------------------------------  //--------------------------------------
623    void TrkLevel2::Set(){
624        if(!Track)Track    = new TClonesArray("TrkTrack");
625        if(!SingletX)SingletX = new TClonesArray("TrkSinglet");
626        if(!SingletY)SingletY = new TClonesArray("TrkSinglet");
627    }
628    //--------------------------------------
629    //
630    //
631    //--------------------------------------
632  void TrkLevel2::Dump(){  void TrkLevel2::Dump(){
633                    
634          //          //
# Line 581  void TrkLevel2::Dump(){ Line 663  void TrkLevel2::Dump(){
663  /**  /**
664   * Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common).   * Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common).
665   */   */
666  void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2){  // void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2){
667    
668    //      //  temporary objects:
669    //     TrkSinglet* t_singlet = new TrkSinglet();
670    //     TrkTrack*   t_track   = new TrkTrack();
671    
672    //      //  **** general variables ****
673    // //    good2 = l2->good2;
674    //     for(Int_t i=0; i<12 ; i++){
675    // //           crc[i] = l2->crc[i];
676    //              good[i] = l2->good[i];
677    //      };
678    //      //  *** TRACKS ***
679    //     if(!Track) Track = new TClonesArray("TrkTrack");
680    //     TClonesArray &t = *Track;
681    //     for(int i=0; i<l2->ntrk; i++){
682    //      t_track->seqno = i;// NBNBNBNB deve sempre essere = i
683    //      t_track->image = l2->image[i]-1;
684    //      //      cout << "track "<<i<<t_track->seqno << t_track->image<<endl;
685    //      t_track->chi2  = l2->chi2_nt[i];
686    //      t_track->nstep = l2->nstep_nt[i];
687    //      for(int it1=0;it1<5;it1++){
688    //          t_track->al[it1] = l2->al_nt[i][it1];
689    //          for(int it2=0;it2<5;it2++)
690    //              t_track->coval[it1][it2] = l2->coval[i][it2][it1];
691    //      };
692    //      for(int ip=0;ip<6;ip++){
693    //          t_track->xgood[ip]  = l2->xgood_nt[i][ip];
694    //          t_track->ygood[ip]  = l2->ygood_nt[i][ip];
695    //          t_track->xm[ip]     = l2->xm_nt[i][ip];
696    //          t_track->ym[ip]     = l2->ym_nt[i][ip];
697    //          t_track->zm[ip]     = l2->zm_nt[i][ip];
698    //          t_track->resx[ip]   = l2->resx_nt[i][ip];
699    //          t_track->resy[ip]   = l2->resy_nt[i][ip];
700    //          t_track->xv[ip]     = l2->xv_nt[i][ip];
701    //          t_track->yv[ip]     = l2->yv_nt[i][ip];
702    //          t_track->zv[ip]     = l2->zv_nt[i][ip];
703    //          t_track->axv[ip]    = l2->axv_nt[i][ip];
704    //          t_track->ayv[ip]    = l2->ayv_nt[i][ip];
705    //          t_track->dedx_x[ip] = l2->dedx_x[i][ip];
706    //          t_track->dedx_y[ip] = l2->dedx_y[i][ip];
707    // //                   t_track->clx[ip] = 0;
708    // //                   t_track->cly[ip] = 0;
709    //      };
710    //      new(t[i]) TrkTrack(*t_track);
711    //      t_track->Clear();
712    //     };
713    // //  *** SINGLETS ***
714    //     if(!SingletX)SingletX = new TClonesArray("TrkSinglet");
715    //     TClonesArray &sx = *SingletX;
716    //     for(int i=0; i<l2->nclsx; i++){
717    //      t_singlet->plane    = l2->planex[i];
718    //      t_singlet->coord[0] = l2->xs[i][0];
719    //      t_singlet->coord[1] = l2->xs[i][1];
720    //      t_singlet->sgnl     = l2->signlxs[i];
721    // //           t_singlet->cls      = 0;
722    //      new(sx[i]) TrkSinglet(*t_singlet);
723    //      t_singlet->Clear();
724    //     }
725    //     if(!SingletY)SingletY = new TClonesArray("TrkSinglet");
726    //     TClonesArray &sy = *SingletY;
727    //     for(int i=0; i<l2->nclsy; i++){
728    //      t_singlet->plane    = l2->planey[i];
729    //      t_singlet->coord[0] = l2->ys[i][0];
730    //      t_singlet->coord[1] = l2->ys[i][1];
731    //      t_singlet->sgnl     = l2->signlys[i];
732    // //           t_singlet->cls      = 0;
733    //      new(sy[i]) TrkSinglet(*t_singlet);
734    //      t_singlet->Clear();
735    //     };
736            
737    //     delete t_track;
738    //     delete t_singlet;
739    // }
740    //--------------------------------------
741    //
742    //
743    //--------------------------------------
744    /**
745     * Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common).
746     * Ref to Level1 data (clusters) is also set. If l1==NULL no references are set.
747     * (NB It make sense to set references only if events are stored in a tree that contains also the Level1 branch)
748     */
749    void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2, TrkLevel1 *l1){
750    
751          //  temporary objects:  //    cout << "void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2, TrkLevel1 *l1)"<<endl;
752        Clear();
753    //  temporary objects:
754      TrkSinglet* t_singlet = new TrkSinglet();      TrkSinglet* t_singlet = new TrkSinglet();
755      TrkTrack*   t_track   = new TrkTrack();      TrkTrack*   t_track   = new TrkTrack();
756    
757          //  **** general variables ****  //  -----------------
758  //    good2 = l2->good2;  //  general variables
759    //  -----------------
760      for(Int_t i=0; i<12 ; i++){      for(Int_t i=0; i<12 ; i++){
761  //              crc[i] = l2->crc[i];          good[i] = l2->good[i];
762                  good[i] = l2->good[i];      };
763          };  //  --------------
764          //  *** TRACKS ***  //  *** TRACKS ***
765    //  --------------
766      if(!Track) Track = new TClonesArray("TrkTrack");      if(!Track) Track = new TClonesArray("TrkTrack");
767      TClonesArray &t = *Track;      TClonesArray &t = *Track;
768        //-----------------------------------------------------
769        if( l1 && !t_track->clx )t_track->clx = new TRefArray(6,0);
770        if( l1 && !t_track->cly )t_track->cly = new TRefArray(6,0);
771        //-----------------------------------------------------
772      for(int i=0; i<l2->ntrk; i++){      for(int i=0; i<l2->ntrk; i++){
773    //      cout <<" TRACK "<<i<<" ------------------ "<<endl;
774          t_track->seqno = i;// NBNBNBNB deve sempre essere = i          t_track->seqno = i;// NBNBNBNB deve sempre essere = i
775          t_track->image = l2->image[i]-1;          t_track->image = l2->image[i]-1;
         //      cout << "track "<<i<<t_track->seqno << t_track->image<<endl;  
776          t_track->chi2  = l2->chi2_nt[i];          t_track->chi2  = l2->chi2_nt[i];
777          t_track->nstep = l2->nstep_nt[i];          t_track->nstep = l2->nstep_nt[i];
778          for(int it1=0;it1<5;it1++){          for(int it1=0;it1<5;it1++){
# Line 608  void TrkLevel2::SetFromLevel2Struct(cTrk Line 781  void TrkLevel2::SetFromLevel2Struct(cTrk
781                  t_track->coval[it1][it2] = l2->coval[i][it2][it1];                  t_track->coval[it1][it2] = l2->coval[i][it2][it1];
782          };          };
783          for(int ip=0;ip<6;ip++){          for(int ip=0;ip<6;ip++){
784              t_track->xgood[ip]  = l2->xgood_nt[i][ip];              t_track->xgood[ip]  = l2->cltrx[i][ip];//l2->xgood_nt[i][ip];
785              t_track->ygood[ip]  = l2->ygood_nt[i][ip];              t_track->ygood[ip]  = l2->cltry[i][ip];//l2->ygood_nt[i][ip];
786              t_track->xm[ip]     = l2->xm_nt[i][ip];              t_track->xm[ip]     = l2->xm_nt[i][ip];
787              t_track->ym[ip]     = l2->ym_nt[i][ip];              t_track->ym[ip]     = l2->ym_nt[i][ip];
788              t_track->zm[ip]     = l2->zm_nt[i][ip];              t_track->zm[ip]     = l2->zm_nt[i][ip];
# Line 622  void TrkLevel2::SetFromLevel2Struct(cTrk Line 795  void TrkLevel2::SetFromLevel2Struct(cTrk
795              t_track->ayv[ip]    = l2->ayv_nt[i][ip];              t_track->ayv[ip]    = l2->ayv_nt[i][ip];
796              t_track->dedx_x[ip] = l2->dedx_x[i][ip];              t_track->dedx_x[ip] = l2->dedx_x[i][ip];
797              t_track->dedx_y[ip] = l2->dedx_y[i][ip];              t_track->dedx_y[ip] = l2->dedx_y[i][ip];
798  //                      t_track->clx[ip] = 0;  //          cout << " ip "<<ip<<" l2->cltrx[i][ip] "<< l2->cltrx[i][ip]<<endl;
799  //                      t_track->cly[ip] = 0;  //          cout << " ip "<<ip<<" l2->cltry[i][ip] "<< l2->cltry[i][ip]<<endl;
800                //-----------------------------------------------------
801                //-----------------------------------------------------
802    //          if(l1 && t_track->xgood[ip])t_track->clx->AddAt(l1->GetCluster(l2->cltrx[i][ip]-1),ip);
803    //          if(l1 && t_track->ygood[ip])t_track->cly->AddAt(l1->GetCluster(l2->cltry[i][ip]-1),ip);
804                if(l2->xgood_nt[i][ip]){
805    //                  cout << " ip "<<ip<<" l2->cltrx[i][ip] "<< l2->cltrx[i][ip]<<endl;
806    //                  cout << " ip "<<ip<<" l2->cltrx[i][ip] "<< l2->cltrx[i][ip]<<" ";
807    //                  if( l1->GetCluster(l2->cltrx[i][ip]-1)->TestBit(l1->GetCluster(l2->cltrx[i][ip]-1)->kIsReferenced) )cout << ">> is referenced ";
808                    if(l1)t_track->clx->AddAt(l1->GetCluster(l2->cltrx[i][ip]-1),ip);
809    //                  cout << " --- "<<l1->GetCluster(l2->cltrx[i][ip]-1)->GetUniqueID()<<endl;
810    //              t_track->xgood[ip] = l2->cltrx[i][ip]; // WORK-AROUND *****
811                }else{
812                    if(l1)t_track->clx->RemoveAt(ip);
813                }
814                if(l2->ygood_nt[i][ip]){
815    //                  cout << " ip "<<ip<<" l2->cltry[i][ip] "<< l2->cltry[i][ip]<<endl;
816    //                  cout << " ip "<<ip<<" l2->cltry[i][ip] "<< l2->cltry[i][ip]<<" ";
817    //                  if( l1->GetCluster(l2->cltry[i][ip]-1)->TestBit(l1->GetCluster(l2->cltry[i][ip]-1)->kIsReferenced) )cout << ">> is referenced ";
818                    if(l1)t_track->cly->AddAt(l1->GetCluster(l2->cltry[i][ip]-1),ip);
819    //                  cout << " --- "<<l1->GetCluster(l2->cltry[i][ip]-1)->GetUniqueID()<<endl;
820    //              t_track->ygood[ip] = l2->cltry[i][ip]; // WORK-AROUND *****
821                }else{
822                    if(l1)t_track->cly->RemoveAt(ip);
823                }
824                //-----------------------------------------------------
825                //-----------------------------------------------------
826          };          };
827          new(t[i]) TrkTrack(*t_track);          new(t[i]) TrkTrack(*t_track);
828          t_track->Clear();          t_track->Clear();
829      };      };
830    //  ----------------
831  //  *** SINGLETS ***  //  *** SINGLETS ***
832    //  ----------------
833      if(!SingletX)SingletX = new TClonesArray("TrkSinglet");      if(!SingletX)SingletX = new TClonesArray("TrkSinglet");
834      TClonesArray &sx = *SingletX;      TClonesArray &sx = *SingletX;
835      for(int i=0; i<l2->nclsx; i++){      for(int i=0; i<l2->nclsx; i++){
# Line 636  void TrkLevel2::SetFromLevel2Struct(cTrk Line 837  void TrkLevel2::SetFromLevel2Struct(cTrk
837          t_singlet->coord[0] = l2->xs[i][0];          t_singlet->coord[0] = l2->xs[i][0];
838          t_singlet->coord[1] = l2->xs[i][1];          t_singlet->coord[1] = l2->xs[i][1];
839          t_singlet->sgnl     = l2->signlxs[i];          t_singlet->sgnl     = l2->signlxs[i];
840  //              t_singlet->cls      = 0;          //-----------------------------------------------------
841            if(l1) t_singlet->cls      = l1->GetCluster(l2->clsx[i]-1);
842            //-----------------------------------------------------
843          new(sx[i]) TrkSinglet(*t_singlet);          new(sx[i]) TrkSinglet(*t_singlet);
844          t_singlet->Clear();          t_singlet->Clear();
845      }      }
# Line 647  void TrkLevel2::SetFromLevel2Struct(cTrk Line 850  void TrkLevel2::SetFromLevel2Struct(cTrk
850          t_singlet->coord[0] = l2->ys[i][0];          t_singlet->coord[0] = l2->ys[i][0];
851          t_singlet->coord[1] = l2->ys[i][1];          t_singlet->coord[1] = l2->ys[i][1];
852          t_singlet->sgnl     = l2->signlys[i];          t_singlet->sgnl     = l2->signlys[i];
853  //              t_singlet->cls      = 0;          //-----------------------------------------------------
854            if(l1) t_singlet->cls      = l1->GetCluster(l2->clsy[i]-1);
855            //-----------------------------------------------------
856          new(sy[i]) TrkSinglet(*t_singlet);          new(sy[i]) TrkSinglet(*t_singlet);
857          t_singlet->Clear();          t_singlet->Clear();
858      };      };
# Line 655  void TrkLevel2::SetFromLevel2Struct(cTrk Line 860  void TrkLevel2::SetFromLevel2Struct(cTrk
860      delete t_track;      delete t_track;
861      delete t_singlet;      delete t_singlet;
862  }  }
 //--------------------------------------  
 //  
 //  
 //--------------------------------------  
 /**  
  * Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common).  
  * Ref to Level1 data (clusters) is also set.  
  */  
 void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2, TrkLevel1 *l1){  
   
 //  temporary objects:  
         TrkSinglet* t_singlet = new TrkSinglet();  
         TrkTrack*   t_track   = new TrkTrack();  
 // general variables  
 //      good2 = l2->good2;  
         for(Int_t i=0; i<12 ; i++){  
 //              crc[i] = l2->crc[i];  
                 good[i] = l2->good[i];  
         };  
 // *** TRACKS ***  
         if(!Track) Track = new TClonesArray("TrkTrack");  
         TClonesArray &t = *Track;  
         if(!t_track->clx)t_track->clx = new TRefArray(6,0);  
         if(!t_track->cly)t_track->cly = new TRefArray(6,0);  
         for(int i=0; i<l2->ntrk; i++){  
                 t_track->seqno = i;// NBNBNBNB deve sempre essere = i  
                 t_track->image = l2->image[i]-1;  
 //              cout << "track "<<i<<t_track->seqno << t_track->image<<endl;  
                 t_track->chi2  = l2->chi2_nt[i];  
                 t_track->nstep = l2->nstep_nt[i];  
                 for(int it1=0;it1<5;it1++){  
                         t_track->al[it1] = l2->al_nt[i][it1];  
                         for(int it2=0;it2<5;it2++)  
                                 t_track->coval[it1][it2] = l2->coval[i][it2][it1];  
                 };  
                 for(int ip=0;ip<6;ip++){  
                         t_track->xgood[ip]  = l2->xgood_nt[i][ip];  
                         t_track->ygood[ip]  = l2->ygood_nt[i][ip];  
                         t_track->xm[ip]     = l2->xm_nt[i][ip];  
                         t_track->ym[ip]     = l2->ym_nt[i][ip];  
                         t_track->zm[ip]     = l2->zm_nt[i][ip];  
                         t_track->resx[ip]   = l2->resx_nt[i][ip];  
                         t_track->resy[ip]   = l2->resy_nt[i][ip];  
                         t_track->xv[ip]     = l2->xv_nt[i][ip];  
                         t_track->yv[ip]     = l2->yv_nt[i][ip];  
                         t_track->zv[ip]     = l2->zv_nt[i][ip];  
                         t_track->axv[ip]    = l2->axv_nt[i][ip];  
                         t_track->ayv[ip]    = l2->ayv_nt[i][ip];  
                         t_track->dedx_x[ip] = l2->dedx_x[i][ip];  
                         t_track->dedx_y[ip] = l2->dedx_y[i][ip];  
 //                      cout << "traccia "<<i<<"  --  "<< ip << " "<< l2->cltrx[i][ip] <<" "<< l2->cltry[i][ip] <<" "<< t_track->xgood[ip] << t_track->ygood[ip]<<endl;  
                         //-----------------------------------------------------  
                         //-----------------------------------------------------  
                         if(t_track->xgood[ip])t_track->clx->AddAt(l1->GetCluster(l2->cltrx[i][ip]-1),ip);  
                         if(t_track->ygood[ip])t_track->cly->AddAt(l1->GetCluster(l2->cltry[i][ip]-1),ip);  
 //                      if(t_track->ygood[ip])cout<<" i "<<i<<" ip "<<ip<<" l2->cltry[i][ip] "<<l2->cltry[i][ip]<< " l1->GetCluster(l2->cltry[i][ip]-1) "<<l1->GetCluster(l2->cltry[i][ip]-1)<<endl;  
 //                      if(t_track->xgood[ip])cout<<" i "<<i<<" ip "<<ip<<" l2->cltrx[i][ip] "<<l2->cltrx[i][ip]<< " l1->GetCluster(l2->cltrx[i][ip]-1) "<<l1->GetCluster(l2->cltrx[i][ip]-1)<<endl;  
                         //-----------------------------------------------------  
                         //-----------------------------------------------------  
                 };  
                 new(t[i]) TrkTrack(*t_track);  
                 t_track->Clear();  
         };  
 // *** SINGLETS ***  
         if(!SingletX)SingletX = new TClonesArray("TrkSinglet");  
         TClonesArray &sx = *SingletX;  
         for(int i=0; i<l2->nclsx; i++){  
                 t_singlet->plane    = l2->planex[i];  
                 t_singlet->coord[0] = l2->xs[i][0];  
                 t_singlet->coord[1] = l2->xs[i][1];  
                 t_singlet->sgnl     = l2->signlxs[i];  
                 //-----------------------------------------------------  
 //              cout << "singolo x "<<i<<"  --  "<<  l2->clsx[i] <<endl;  
                 t_singlet->cls      = l1->GetCluster(l2->clsx[i]-1);  
 //              cout<<" i "<<i<<" l2->clsx[i] "<<l2->clsx[i]<< " l1->GetCluster(l2->clsx[i]-1) "<<l1->GetCluster(l2->clsx[i]-1)<<endl;            
                 //-----------------------------------------------------  
                 new(sx[i]) TrkSinglet(*t_singlet);  
                 t_singlet->Clear();  
         }  
         if(!SingletY)SingletY = new TClonesArray("TrkSinglet");  
         TClonesArray &sy = *SingletY;  
         for(int i=0; i<l2->nclsy; i++){  
                 t_singlet->plane    = l2->planey[i];  
                 t_singlet->coord[0] = l2->ys[i][0];  
                 t_singlet->coord[1] = l2->ys[i][1];  
                 t_singlet->sgnl     = l2->signlys[i];  
                 //-----------------------------------------------------  
 //              cout << "singolo y "<<i<<"  --  "<<  l2->clsy[i] <<endl;  
                 t_singlet->cls      = l1->GetCluster(l2->clsy[i]-1);  
 //              cout<<" i "<<i<<" l2->clsy[i] "<<l2->clsy[i]<< " l1->GetCluster(l2->clsy[i]-1) "<<l1->GetCluster(l2->clsy[i]-1)<<endl;            
                 //-----------------------------------------------------  
                 new(sy[i]) TrkSinglet(*t_singlet);  
                 t_singlet->Clear();  
         };  
           
         delete t_track;  
         delete t_singlet;  
 }  
863  /**  /**
864   * Fills a struct cTrkLevel2 with values from a TrkLevel2 object (to put data into a F77 common).   * Fills a struct cTrkLevel2 with values from a TrkLevel2 object (to put data into a F77 common).
865   */   */
# Line 779  void TrkLevel2::GetLevel2Struct(cTrkLeve Line 886  void TrkLevel2::GetLevel2Struct(cTrkLeve
886                      l2->coval[i][it2][it1] = ((TrkTrack *)Track->At(i))->coval[it1][it2];                      l2->coval[i][it2][it1] = ((TrkTrack *)Track->At(i))->coval[it1][it2];
887              };              };
888              for(int ip=0;ip<6;ip++){              for(int ip=0;ip<6;ip++){
889                  l2->xgood_nt[i][ip] = ((TrkTrack *)Track->At(i))->xgood[ip];                  l2->xgood_nt[i][ip] = ((TrkTrack *)Track->At(i))->XGood(ip);
890                  l2->ygood_nt[i][ip] = ((TrkTrack *)Track->At(i))->ygood[ip];                  l2->ygood_nt[i][ip] = ((TrkTrack *)Track->At(i))->YGood(ip);
891                  l2->xm_nt[i][ip]    = ((TrkTrack *)Track->At(i))->xm[ip];                  l2->xm_nt[i][ip]    = ((TrkTrack *)Track->At(i))->xm[ip];
892                  l2->ym_nt[i][ip]    = ((TrkTrack *)Track->At(i))->ym[ip];                  l2->ym_nt[i][ip]    = ((TrkTrack *)Track->At(i))->ym[ip];
893                  l2->zm_nt[i][ip]    = ((TrkTrack *)Track->At(i))->zm[ip];                  l2->zm_nt[i][ip]    = ((TrkTrack *)Track->At(i))->zm[ip];
# Line 867  TRefArray *TrkLevel2::GetTracks_NFitSort Line 974  TRefArray *TrkLevel2::GetTracks_NFitSort
974                    
975      int indo=0;      int indo=0;
976      int indi=0;      int indi=0;
977      while(N != 0){      while(N > 0){
978    //    while(N != 0){
979          int nfit =0;          int nfit =0;
980          float chi2ref = numeric_limits<float>::max();          float chi2ref = numeric_limits<float>::max();
981                                    
# Line 878  TRefArray *TrkLevel2::GetTracks_NFitSort Line 986  TRefArray *TrkLevel2::GetTracks_NFitSort
986              }              }
987          }          }
988          //second loop to search minimum chi2 among selected          //second loop to search minimum chi2 among selected
989          for(int i=0; i<this->ntrk(); i++){          for(int i=0; i<ntrk(); i++){
990              Float_t chi2 = ((TrkTrack *)t[i])->chi2;              Float_t chi2 = ((TrkTrack *)t[i])->chi2;
991              if(chi2 < 0) chi2 = chi2*1000;              if(chi2 < 0) chi2 = -chi2*1000;
992              if(    chi2 < chi2ref              if(    chi2 < chi2ref
993                     && ((TrkTrack *)t[i])->GetNtot() == nfit                     && ((TrkTrack *)t[i])->GetNtot() == nfit
994                     && m[i]==1){                     && m[i]==1){
# Line 892  TRefArray *TrkLevel2::GetTracks_NFitSort Line 1000  TRefArray *TrkLevel2::GetTracks_NFitSort
1000              m[((TrkTrack *)t[indi])->image] = 0;              m[((TrkTrack *)t[indi])->image] = 0;
1001              N--;              N--;
1002                    
1003              //      cout << "i** "<< ((TrkTrack *)t[indi])->image << " " << nfiti <<" "<<chi2i<<endl;  //          cout << "i** "<< ((TrkTrack *)t[indi])->image << " " << nfiti <<" "<<chi2i<<endl;
1004          };          };
1005          sorted->Add( (TrkTrack*)t[indi] );                sorted->Add( (TrkTrack*)t[indi] );      
1006                                    
1007          m[indi] = 0;          m[indi] = 0;
1008  //              cout << "SORTED "<< indo << " "<< indi << " "<< N << endl;  //      cout << "SORTED "<< indo << " "<< indi << " "<< N << " "<<((TrkTrack *)t[indi])->image<<" "<<chi2ref<<endl;
1009          N--;              N--;    
1010          indo++;          indo++;
1011      }      }
1012      m.clear();      m.clear();
1013  //      cout << "GetTracks_NFitSorted(it): Done"<< endl;  //    cout << "GetTracks_NFitSorted(it): Done"<< endl;
1014    
1015      return sorted;      return sorted;
1016  //    return PhysicalTrack;  //    return PhysicalTrack;
# Line 1052  TrkTrack *TrkLevel2::GetTrackImage(int i Line 1160  TrkTrack *TrkLevel2::GetTrackImage(int i
1160   */   */
1161  void TrkLevel2::LoadField(TString path){  void TrkLevel2::LoadField(TString path){
1162  //  //
1163      strcpy(path_.path,path.Data());  //     strcpy(path_.path,path.Data());
1164      path_.pathlen = path.Length();  //     path_.pathlen = path.Length();
1165      path_.error   = 0;  //     path_.error   = 0;
1166      readb_();  //     readb_();
1167    
1168        TrkParams::Set(path,1);
1169        TrkParams::Load(1);
1170    
1171  //  //
1172  };  };
1173    /**
1174     * Get BY (kGauss)
1175     * @param v (x,y,z) coordinates in cm
1176     */
1177    float TrkLevel2::GetBX(float* v){
1178        float b[3];
1179        gufld_(v,b);
1180        return b[0]/10.;
1181    }
1182    /**
1183     * Get BY (kGauss)
1184     * @param v (x,y,z) coordinates in cm
1185     */
1186    float TrkLevel2::GetBY(float* v){
1187        float b[3];
1188        gufld_(v,b);
1189        return b[1]/10.;
1190    }
1191    /**
1192     * Get BY (kGauss)
1193     * @param v (x,y,z) coordinates in cm
1194     */
1195    float TrkLevel2::GetBZ(float* v){
1196        float b[3];
1197        gufld_(v,b);
1198        return b[2]/10.;
1199    }
1200  //--------------------------------------  //--------------------------------------
1201  //  //
1202  //  //
# Line 1235  int Trajectory::DoTrack2(float* al){ Line 1374  int Trajectory::DoTrack2(float* al){
1374      for (int i=0; i<5; i++)      dal[i]  = (double)al[i];      for (int i=0; i<5; i++)      dal[i]  = (double)al[i];
1375      for (int i=0; i<npoint; i++) dzin[i] = (double)z[i];      for (int i=0; i<npoint; i++) dzin[i] = (double)z[i];
1376    
1377        TrkParams::Load(1);
1378        if( !TrkParams::IsLoaded(1) ){
1379            cout << "int Trajectory::DoTrack2(float* al) --- ERROR --- m.field not loaded"<<endl;
1380            return 0;
1381        }
1382      dotrack2_(&(npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail);      dotrack2_(&(npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail);
1383            
1384      for (int i=0; i<npoint; i++){      for (int i=0; i<npoint; i++){
Legend:
Removed from v.1.21  
changed lines
  Added in v.1.31
  ViewVC Help
Powered by ViewVC 1.1.23