/[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.1 by mocchiut, Fri May 19 13:15:54 2006 UTC revision 1.2 by pam-fi, Thu Jun 1 09:03:09 2006 UTC
# Line 10  using namespace std; Line 10  using namespace std;
10  //......................................  //......................................
11  extern "C" {      extern "C" {    
12      void dotrack_(int*, double*, double*, double*, double*, int*);      void dotrack_(int*, double*, double*, double*, double*, int*);
13        void dotrack2_(int*, double*, double*, double*, double*,double*, double*, double*,int*);
14      int  readb_(const char*);      int  readb_(const char*);
15  }  }
16  //--------------------------------------  //--------------------------------------
# Line 110  int TrkTrack::DoTrack(Trajectory* t){ Line 111  int TrkTrack::DoTrack(Trajectory* t){
111  //  //
112  //  //
113  //--------------------------------------  //--------------------------------------
114    /**
115     * Evaluates the trajectory in the apparatus associated to the track.
116     * It integrates the equations of motion in the magnetic field. The magnetic field should be previously loaded ( by calling  TrkLevel2::LoadField() ), otherwise an error message is returned.  
117     * @param t pointer to an object of the class Trajectory,
118     * which z coordinates should be previously initialized by calling the proper constructor ( Trajectory::Trajectory(int n, float* zin) ).
119     * @return error flag.
120     */
121    int TrkTrack::DoTrack2(Trajectory* t){
122    
123        double *dxout   = new double[t->npoint];
124        double *dyout   = new double[t->npoint];
125        double *dthxout = new double[t->npoint];
126        double *dthyout = new double[t->npoint];
127        double *dtlout  = new double[t->npoint];
128        double *dzin    = new double[t->npoint];
129        double dal[5];
130    
131        int ifail = 0;
132    
133        for (int i=0; i<5; i++)         dal[i]  = (double)al[i];
134        for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i];
135    
136        dotrack2_(&(t->npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail);
137        
138        for (int i=0; i<t->npoint; i++){
139            t->x[i]   = (float)*dxout++;
140            t->y[i]   = (float)*dyout++;
141            t->thx[i] = (float)*dthxout++;
142            t->thy[i] = (float)*dthyout++;
143            t->tl[i]  = (float)*dtlout++;
144        }
145    
146    //    delete [] dxout;
147    //    delete [] dyout;
148    //    delete [] dzin;
149    
150        return ifail;
151    };
152    //--------------------------------------
153    //
154    //
155    //--------------------------------------
156  //float TrkTrack::BdL(){  //float TrkTrack::BdL(){
157  //};  //};
158  //--------------------------------------  //--------------------------------------
# Line 408  void TrkLevel2::LoadField(TString s){ Line 451  void TrkLevel2::LoadField(TString s){
451  //  //
452  //--------------------------------------  //--------------------------------------
453  /**  /**
454     * Trajectory default constructor.
455     * (By default is created with z-coordinates inside the tracking volume)
456      */
457    Trajectory::Trajectory(){
458        npoint = 10;
459        x = new float[npoint];
460        y = new float[npoint];
461        z = new float[npoint];
462        thx = new float[npoint];
463        thy = new float[npoint];
464        tl = new float[npoint];
465        float dz = ((ZTRKUP)-(ZTRKDW))/(npoint-1);
466        for(int i=0; i<npoint; i++){
467            x[i] = 0;
468            y[i] = 0;
469            z[i] = (ZTRKUP) - i*dz;
470            thx[i] = 0;
471            thy[i] = 0;
472            tl[i] = 0;
473        }
474    }
475    //--------------------------------------
476    //
477    //
478    //--------------------------------------
479    /**
480   * Trajectory constructor.   * Trajectory constructor.
481     * (By default is created with z-coordinates inside the tracking volume)
482   * \param n Number of points   * \param n Number of points
483   */   */
484  Trajectory::Trajectory(int n){  Trajectory::Trajectory(int n){
485        if(n<=0){
486            cout << "NB! Trajectory must have at least 1 point >>> created with 10 points" << endl;
487            n=10;
488        }
489      npoint = n;      npoint = n;
490      x = new float[npoint];      x = new float[npoint];
491      y = new float[npoint];      y = new float[npoint];
492      z = new float[npoint];      z = new float[npoint];
493        thx = new float[npoint];
494        thy = new float[npoint];
495        tl = new float[npoint];
496        float dz = ((ZTRKUP)-(ZTRKDW))/(npoint-1);
497      for(int i=0; i<npoint; i++){      for(int i=0; i<npoint; i++){
498          x[i] = 0;          x[i] = 0;
499          y[i] = 0;          y[i] = 0;
500          z[i] = 0;          z[i] = (ZTRKUP) - i*dz;
501            thx[i] = 0;
502            thy[i] = 0;
503            tl[i] = 0;
504      }      }
505  }  }
506  //--------------------------------------  //--------------------------------------
# Line 432  Trajectory::Trajectory(int n){ Line 513  Trajectory::Trajectory(int n){
513   * \param pz Pointer to float array, defining z coordinates   * \param pz Pointer to float array, defining z coordinates
514   */   */
515  Trajectory::Trajectory(int n, float* zin){  Trajectory::Trajectory(int n, float* zin){
516      npoint = n;      npoint = 10;
517        if(n>0)npoint = n;
518      x = new float[npoint];      x = new float[npoint];
519      y = new float[npoint];      y = new float[npoint];
520      z = new float[npoint];      z = new float[npoint];
521      for(int i=0; i<npoint; i++){      thx = new float[npoint];
522        thy = new float[npoint];
523        tl = new float[npoint];
524        int i=0;
525        do{
526          x[i] = 0;          x[i] = 0;
527          y[i] = 0;          y[i] = 0;
528          z[i] = zin[i];          z[i] = zin[i];
529      }          thx[i] = 0;
530            thy[i] = 0;
531            tl[i] = 0;
532            i++;            
533        }while(zin[i-1] > zin[i] && i < npoint);
534        npoint=i;
535        if(npoint != n)cout << "NB! Trajectory created with "<<npoint<<" points"<<endl;
536  }  }
537  //--------------------------------------  //--------------------------------------
538  //  //
# Line 452  Trajectory::Trajectory(int n, float* zin Line 544  Trajectory::Trajectory(int n, float* zin
544  void Trajectory::Dump(){  void Trajectory::Dump(){
545      cout <<endl<< "Trajectory ========== "<<endl;      cout <<endl<< "Trajectory ========== "<<endl;
546      for (int i=0; i<npoint; i++){      for (int i=0; i<npoint; i++){
547          cout << i <<" >> " << x[i] <<" "<< y[i] <<" "<< z[i] << endl;;          cout << i <<" >> " << x[i] <<" "<< y[i] <<" "<< z[i] ;
548            cout <<" -- " << thx[i] <<" "<< thy[i] ;
549            cout <<" -- " << tl[i] << endl;
550      };      };
551  }  }
552    //--------------------------------------
553    //
554    //
555    //--------------------------------------
556    /**
557     * Get trajectory length between two points
558     * @param ifirst first point (default 0)
559     * @param ilast last point (default npoint)
560     */
561    float Trajectory::GetLength(int ifirst, int ilast){
562        if( ifirst<0    ) ifirst = 0;
563        if( ilast>=npoint) ilast  = npoint-1;
564        float l=0;
565        for(int i=ifirst;i<=ilast;i++){
566            l=l+tl[i];
567        };
568        if(z[ilast] > ZINI)l=l-tl[ilast];
569        if(z[ifirst] < ZINI)   l=l-tl[ifirst];
570        
571        return l;
572    
573    }
574  ClassImp(TrkLevel2);  ClassImp(TrkLevel2);
575  ClassImp(TrkSinglet);  ClassImp(TrkSinglet);
576  ClassImp(TrkTrack);  ClassImp(TrkTrack);
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.2
  ViewVC Help
Powered by ViewVC 1.1.23