/[PAMELA software]/PamelaLevel2/doc/examples/example2.C
ViewVC logotype

Diff of /PamelaLevel2/doc/examples/example2.C

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

revision 1.1.1.1 by pam-fi, Fri Jun 16 16:43:55 2006 UTC revision 1.3 by pam-fi, Mon Jan 15 11:51:39 2007 UTC
# Line 4  Line 4 
4  //      The tracking method evaluates, besides the track intersection coordinates at given z-coordinates,  //      The tracking method evaluates, besides the track intersection coordinates at given z-coordinates,
5  //      the track length (total and between two given points along the trajectory) and the projected angles.  //      the track length (total and between two given points along the trajectory) and the projected angles.
6  //  //
7    //  --- Modified on January 2007 ---
8    //
9    //      see the comment at the beginning of example1.C      
10    //
11  example2(TString file){  example2(TString file){
12                    
13          gROOT->Reset();      gROOT->Reset();
14          //      //
15          // create some histograms      // create some histograms
16          //      //
17          TH1F* s1x = new TH1F("s1x","Track x impact position on S1",100,-30.,30.);      TH1F* s1x = new TH1F("s1x","Track x impact position on S1",100,-30.,30.);
18          TH1F* s1y = new TH1F("s1y","Track y impact position on S1",100,-30.,30.);      TH1F* s1y = new TH1F("s1y","Track y impact position on S1",100,-30.,30.);
19          TH1F* s2x = new TH1F("s2x","Track x impact position on S2",100,-30.,30.);      TH1F* s2x = new TH1F("s2x","Track x impact position on S2",100,-30.,30.);
20          TH1F* s2y = new TH1F("s2y","Track y impact position on S2",100,-30.,30.);      TH1F* s2y = new TH1F("s2y","Track y impact position on S2",100,-30.,30.);
21          TH1F* s3x = new TH1F("s3x","Track x impact position on S3",100,-30.,30.);      TH1F* s3x = new TH1F("s3x","Track x impact position on S3",100,-30.,30.);
22          TH1F* s3y = new TH1F("s3y","Track y impact position on S3",100,-30.,30.);      TH1F* s3y = new TH1F("s3y","Track y impact position on S3",100,-30.,30.);
23                    
24          TH1F* s1tx = new TH1F("s1tx","Track x projected angle on S1",50,-20.,20.);      TH1F* s1tx = new TH1F("s1tx","Track x projected angle on S1",50,-20.,20.);
25          TH1F* s1ty = new TH1F("s1ty","Track y projected angle on S1",50,-20.,20.);      TH1F* s1ty = new TH1F("s1ty","Track y projected angle on S1",50,-20.,20.);
26          TH1F* s2tx = new TH1F("s2tx","Track x projected angle on S2",50,-20.,20.);      TH1F* s2tx = new TH1F("s2tx","Track x projected angle on S2",50,-20.,20.);
27          TH1F* s2ty = new TH1F("s2ty","Track y projected angle on S2",50,-20.,20.);      TH1F* s2ty = new TH1F("s2ty","Track y projected angle on S2",50,-20.,20.);
28          TH1F* s3tx = new TH1F("s3tx","Track x projected angle on S3",50,-20.,20.);      TH1F* s3tx = new TH1F("s3tx","Track x projected angle on S3",50,-20.,20.);
29          TH1F* s3ty = new TH1F("s3ty","Track y projected angle on S3",50,-20.,20.);            TH1F* s3ty = new TH1F("s3ty","Track y projected angle on S3",50,-20.,20.);  
30                    
31          TH1F* trl   = new TH1F("trl","Track length",100,0.,100.);      TH1F* trl   = new TH1F("trl","Track length",100,0.,100.);
32          TH1F* trl12 = new TH1F("trl12","Track length S1-S2",100,0.,100.);      TH1F* trl12 = new TH1F("trl12","Track length S1-S2",100,0.,100.);
33          TH1F* trl23 = new TH1F("trl23","Track length S2-S3",100,0.,100.);      TH1F* trl23 = new TH1F("trl23","Track length S2-S3",100,0.,100.);
34                    
35                    
36          //      //
37          PamLevel2*  pam_event  = new PamLevel2();       // << create pamela event      PamLevel2*  pam_event  = new PamLevel2();       // << create pamela event
38          //      //
39          TFile f(file);      TFile f(file);
40          //      //
41          TTree *T = pam_event->LoadPamTrees(&f);         // << load Pamela trees from file f      TTree *T = pam_event->GetPamTree(&f);         // << get Pamela trees from file f
42          Int_t nevent = T->GetEntries();      Int_t nevent = T->GetEntries();    
43          // ********************      // ********************
44          // load magnetic field      // load magnetic field
45          // ********************      // ********************
46          pam_event->LoadField("./magnetic-field/");      //      pam_event->GetTrkLevel2()->LoadField("./magnetic-field/");  //
47          // ********************      // ********************
48          // initialize some trajectories      // initialize some trajectories
49          // ********************      // ********************
50          // << default trajectory is created with 10 points between the upper and lower tracker planes      // << default trajectory is created with 10 points between the upper and lower tracker planes
51          Trajectory *tr1 = new Trajectory() ;                                   // << create a default trajectory      Trajectory *tr1 = new Trajectory() ;                                   // << create a default trajectory
52          // << the number of points can be set by the user      // << the number of points can be set by the user
53          Trajectory *tr2 = new Trajectory(100) ;                                // << create a trajectory with 100 points      Trajectory *tr2 = new Trajectory(100) ;                                // << create a trajectory with 100 points
54          // << and also the z-coordinates can be set by the user.      // << and also the z-coordinates can be set by the user.
55          // << for example if we want to evaluate the track intersection points in the TOF planes      // << for example if we want to evaluate the track intersection points in the TOF planes
56          // << we can define the following trajectory      // << we can define the following trajectory
57          Int_t nz = 6; Float_t zin[6];                                          // << define TOF z-coordinates      Int_t nz = 6; Float_t zin[6];                                          // << define TOF z-coordinates
58          for(Int_t ip=0; ip<nz; ip++)      for(Int_t ip=0; ip<nz; ip++)
59                 zin[ip] = pam_event->GetZTOF(pam_event->GetToFPlaneID(ip));     // << read ToF plane z-coordinates          zin[ip] = pam_event->GetToFLevel2()->GetZTOF(pam_event->GetToFLevel2()->GetToFPlaneID(ip));     // << read ToF plane z-coordinates
60          Trajectory *tr = new Trajectory(nz,zin);                               // << create a trajectory in the apparatus      Trajectory *tr = new Trajectory(nz,zin);                           // << create a trajectory in the apparatus
61          //      //
62          cout << endl<< " Start loop over events   ";      cout << endl<< " Start loop over events   ";
63          for (Int_t i=0; i<nevent;i++){      for (Int_t i=0; i<nevent;i++){
64                  //          //
65                  T->GetEntry(i);          pam_event->Clear();
66                  //  
67                  if(pam_event->GetNTracks()==1){                 // << select events with only one track          T->GetEntry(i);
68                      //          //
69                      TrkTrack *track = pam_event->GetTrack(0);   // << retrieve the track          if(pam_event->GetTrkLevel2()->GetNTracks()==1){                 // << select events with only one track
70                          //              //
71                          // << perform some track selection              PamTrack *track = pam_event->GetTrack(0);   // << retrieve the track
72                          //              //
73                          if(                                                  // << perform some track selection
74                                  track->chi2 > 0    &&              //
75                                  track->chi2 < 100  &&              if(                                
76                                  track->GetNX() >= 4  &&                  track->GetTrkTrack()->chi2 > 0    &&
77                                  track->GetNY() >= 3  &&                  track->GetTrkTrack()->chi2 < 100  &&
78                                  true                  track->GetTrkTrack()->GetNX() >= 4  &&
79                                  ){                        track->GetTrkTrack()->GetNY() >= 3  &&
80                    true
81                    ){      
82                                                                    
83                              cout << endl<< "***** First trajectory"<< endl;                  cout << endl<< "***** First trajectory"<< endl;
84                              track->DoTrack2(tr1);       // << calculate the first trajectory in magnetic field                  track->GetTrkTrack()->DoTrack2(tr1);    // << calculate the first trajectory in magnetic field
85  //                          tr1->Dump();                // dump the trajectory  //                          tr1->Dump();                // dump the trajectory
86                              cout << "Length: "<< tr1->GetLength()<<endl;  // << get the track length                  cout << "Length: "<< tr1->GetLength()<<endl;  // << get the track length
87                              //                  //
88                              cout << endl<< "***** Second trajectory"<< endl;                  cout << endl<< "***** Second trajectory"<< endl;
89                              track->DoTrack2(tr2);       // << calculate trajectory in magnetic field                  track->GetTrkTrack()->DoTrack2(tr2);    // << calculate trajectory in magnetic field
90                              cout << "Length: "<< tr2->GetLength()<<endl; // << get the track length                  cout << "Length: "<< tr2->GetLength()<<endl; // << get the track length
91                              //                  //
92                              cout << endl<< "***** Third trajectory"<< endl;                  cout << endl<< "***** Third trajectory"<< endl;
93                              track->DoTrack2(tr);        // << calculate trajectory in magnetic field                                              track->GetTrkTrack()->DoTrack2(tr);     // << calculate trajectory in magnetic field                            
94  //                            tr->Dump();  //                            tr->Dump();
95                              cout << "Length (S11-S32): "<< tr->GetLength()<<endl;                  cout << "Length (S11-S32): "<< tr->GetLength()<<endl;
96                              // << The length can be evaluated also bewteen two planes set by the user                  // << The length can be evaluated also bewteen two planes set by the user
97                              cout << "Length between S11-S21: "<< tr->GetLength(0,2)<< endl;                  cout << "Length between S11-S21: "<< tr->GetLength(0,2)<< endl;
98                              cout << "Length between S21-S31: "<< tr->GetLength(2,4)<< endl;                  cout << "Length between S21-S31: "<< tr->GetLength(2,4)<< endl;
99                              cout << "Length between S11-S31: "<< tr->GetLength(0,4)<< endl;                  cout << "Length between S11-S31: "<< tr->GetLength(0,4)<< endl;
100                              // now fills some histos:                  // now fills some histos:
101                              // store calculated coordinates                  // store calculated coordinates
102                              s1x->Fill( tr->x[0] );                        s1x->Fill( tr->x[0] );  
103                              s1y->Fill( tr->y[0] );                  s1y->Fill( tr->y[0] );
104                              s2x->Fill( tr->x[1] );                  s2x->Fill( tr->x[1] );
105                              s2y->Fill( tr->y[1] );                  s2y->Fill( tr->y[1] );
106                              s3x->Fill( tr->x[2] );                  s3x->Fill( tr->x[2] );
107                              s3y->Fill( tr->y[2] );                  s3y->Fill( tr->y[2] );
108                              // store calculated projected angles                  // store calculated projected angles
109                              s1tx->Fill( tr->thx[0] );                    s1tx->Fill( tr->thx[0] );      
110                              s1ty->Fill( tr->thy[0] );                  s1ty->Fill( tr->thy[0] );
111                              s2tx->Fill( tr->thx[1] );                  s2tx->Fill( tr->thx[1] );
112                              s2ty->Fill( tr->thy[1] );                  s2ty->Fill( tr->thy[1] );
113                              s3tx->Fill( tr->thx[2] );                  s3tx->Fill( tr->thx[2] );
114                              s3ty->Fill( tr->thy[2] );                  s3ty->Fill( tr->thy[2] );
115                                                            
116                              // store calculated track lengths                  // store calculated track lengths
117                              trl->Fill( tr->GetLength() );                  trl->Fill( tr->GetLength() );
118                              trl12->Fill( tr->GetLength(0,1) );                  trl12->Fill( tr->GetLength(0,1) );
119                              trl23->Fill( tr->GetLength(1,2) );                  trl23->Fill( tr->GetLength(1,2) );
120                              };              };
                 };  
121          };          };
122          cout << endl << " Done "<< endl<<endl;      };
123          //      cout << endl << " Done "<< endl<<endl;
124          // close file      //
125          //      // close file
126          f.Close();            //
127          //      f.Close();  
128          // plot histos      //
129          //      // plot histos
130          }        //
131    }  
Legend:
Removed from v.1.1.1.1  
changed lines
  Added in v.1.3
  ViewVC Help
Powered by ViewVC 1.1.23