ROCsoft: beamtime/cosy-jan12/go4/TTrackingProc.cxx (r4864/r3160) Source File

Go to the documentation of this file.
00001 #include "TTrackingProc.h"
00002 #include "TGo4Version.h"
00003 #if __GO4BUILDVERSION__ > 40502
00004 #include "go4iostream.h"
00005 #else
00006 #include "Riostream.h"
00007 #endif
00008 #include "TTimeStamp.h"
00009 #include "TSystem.h"
00010 #include "TROOT.h"
00011 #include "TGo4Log.h"
00012 
00013 #include "TGo4UserException.h"
00014 
00015 #include "STS/TSTSGeometry.h"
00016 #include "TSpacepoint.h"
00017 #include "STS/TSTSParam.h"
00018 
00019 
00020 TTrackingProc::TTrackingProc(const char* name) :
00021    TCBMBeamtimeProc(name),
00022    fBeamMonitorEvent(0),
00023    fFirstHodoEvent(0),
00024    fSecondHodoEvent(0),
00025    fGEMEvent(0), 
00026    fSTSEvent(0),
00027    fOutputEvent(0)
00028 {
00029    cout << "**** TTrackingProc: Create instance ****" << name << endl;
00030    TString obname, obtitle;
00031 
00032    fPar = (TTrackingParam*) MakeParameter("TrackingPar", "TTrackingParam");
00033 
00034 // fKalman = new CbmKF();
00035 
00036   fPosXHod1 = MakeTH1('I', "Tracking/XPosHod1","XPosHod1_title",  1000,-50.,50.);
00037   fPosYHod1 = MakeTH1('I', "Tracking/YPosHod1","YPosHod1_title",  1000,-50.,50.);
00038   fPosXYHod1 = MakeTH2('I', "Tracking/XYPosHod1","XYPosHod1_title",  100,-50.,50.,  100,-50.,50.);
00039   fPosXYHod1Sts = MakeTH2('I', "Tracking/XYPosHod1Sts","XYPosHod1Sts_title",  100,-50.,50.,  100,-50.,50.);
00040   fPosXHod2 = MakeTH1('I', "Tracking/XPosHod2","XPosHod2_title",  1000,-50.,50.);
00041   fPosYHod2 = MakeTH1('I', "Tracking/YPosHod2","YPosHod2_title",  1000,-50.,50.);
00042   fPosXYHod2 = MakeTH2('I', "Tracking/XYPosHod2","XYPosHod2_title",  100,-50.,50.,  100,-50.,50.);
00043   fPosXYHod2Sts = MakeTH2('I', "Tracking/XYPosHod2Sts","XYPosHod2Sts_title",  100,-50.,50.,  100,-50.,50.);
00044 
00045   fResXHod1 = MakeTH1('I', "Tracking/XResHod1","XResHod1_title",  400,-20.,20.);
00046   fResYHod1 = MakeTH1('I', "Tracking/YResHod1","YResHod1_title",  400,-20.,20.);
00047   fResXHod2 = MakeTH1('I', "Tracking/XResHod2","XResHod2_title",  400,-20.,20.);
00048   fResYHod2 = MakeTH1('I', "Tracking/YResHod2","YResHod2_title",  400,-20.,20.);
00049 
00050   fCorrYHod = MakeTH2('I', "Tracking/CorrYHod","CorrYHod_title",  100,-50.,50.,  100,-50.,50.); 
00051   fCorrXHod = MakeTH2('I', "Tracking/CorrXHod","CorrXHod_title",  100,-50.,50.,  100,-50.,50.); 
00052   fCorrYHod1Gem1 = MakeTH2('I', "Tracking/CorrYHod1Gem1","CorrYHod1Gem1_title",  100,-50.,50.,  100,-50.,50.); 
00053   fCorrXHod1Gem1 = MakeTH2('I', "Tracking/CorrXHod1Gem1","CorrXHod1Gem1_title",  100,-50.,50.,  100,-50.,50.); 
00054   fCorrYHod2Gem1 = MakeTH2('I', "Tracking/CorrYHod2Gem1","CorrYHod2Gem1_title",  100,-50.,50.,  100,-50.,50.); 
00055   fCorrXHod2Gem1 = MakeTH2('I', "Tracking/CorrXHod2Gem1","CorrXHod2Gem1_title",  100,-50.,50.,  100,-50.,50.); 
00056 
00057   TSTSParam * stsParam = dynamic_cast< TSTSParam * >( GetParameter( "STSPar" ) );
00058   if( ! stsParam ) {
00059       printf("Error in TTrackingProc::TTrackingProc: StsPar parameter not found\n");
00060       return;
00061   }
00062   UInt_t detArrSize = stsParam->fSTSTopology->GetDetArrSize();
00063   
00064   fPosXSts.resize( detArrSize );
00065   fPosYSts.resize( detArrSize );
00066   fPosXYSts.resize( detArrSize );
00067   fCorrXXSts.resize( detArrSize );
00068   fCorrYYSts.resize( detArrSize );
00069   fCorrXXH1Sts.resize( detArrSize );
00070   fCorrYYH1Sts.resize( detArrSize );
00071   fCorrXXGem1Sts.resize( detArrSize );
00072   fCorrYYGem1Sts.resize( detArrSize );
00073   fNHitsSts.resize( detArrSize );
00074   fPosXSts_vs_time.resize( detArrSize );
00075   fPosYSts_vs_time.resize( detArrSize );
00076 
00077   for( UInt_t iDet = 0; iDet < detArrSize; iDet++ ) {
00078 //        fPosXSts.at( iDet ) = MakeTH1('I', Form( "Tracking/XPosSts%d", iDet), Form("XPosSts%d_title", iDet),  140,-7.,7.);
00079 //        fPosYSts.at( iDet ) = MakeTH1('I', Form( "Tracking/YPosSts%d", iDet), Form("XPosSts%d_title", iDet),  140,-7.,7.);
00080 //        fPosXYSts.at( iDet ) = MakeTH2('I', Form( "Tracking/XYPosSts%d", iDet), Form("XPosSts%d_title", iDet),  140,-7.,7.,  140,-7.,7. );
00081           fPosXSts.at( iDet ) = MakeTH1('I', Form( "Tracking/XPosSts%d", iDet), Form("XPosSts%d_title", iDet),  10000,-10.,10.);
00082           fPosYSts.at( iDet ) = MakeTH1('I', Form( "Tracking/YPosSts%d", iDet), Form("XPosSts%d_title", iDet),  10000,-10.,10.);
00083           fPosXYSts.at( iDet ) = MakeTH2('I', Form( "Tracking/XYPosSts%d", iDet), Form("XPosSts%d_title", iDet),  100,-10.,10.,  100,-10.,10. );
00084           fCorrXXSts.at( iDet ) = MakeTH2('I', Form( "Tracking/XXCorrSts%d", iDet), Form("XXCorrSts%d_title", iDet),  100,-10.,10.,  100,-10.,10. );
00085           fCorrYYSts.at( iDet ) = MakeTH2('I', Form( "Tracking/YYCorrSts%d", iDet), Form("YYCorrSts%d_title", iDet),  100,-10.,10.,  100,-10.,10. );
00086           fCorrXXH1Sts.at( iDet ) = MakeTH2('I', Form( "Tracking/XXCorrH1Sts%d", iDet), Form("XXCorrH1Sts%d_title", iDet),  250,-20.,30.,  100,-10.,10. );
00087           fCorrYYH1Sts.at( iDet ) = MakeTH2('I', Form( "Tracking/YYCorrH1Sts%d", iDet), Form("YYCorrH1Sts%d_title", iDet),  100,-10.,10.,  100,-10.,10. );
00088           fCorrXXGem1Sts.at( iDet ) = MakeTH2('I', Form( "Tracking/XXCorrGem1Sts%d", iDet), Form("XXCorrGem1Sts%d_title", iDet),  100,-50.,50.,  100,-50.,50. );
00089           fCorrYYGem1Sts.at( iDet ) = MakeTH2('I', Form( "Tracking/YYCorrGem1Sts%d", iDet), Form("YYCorrGem1Sts%d_title", iDet),  100,-50.,50.,  100,-50.,50. );
00090           fNHitsSts.at( iDet ) = MakeTH1('I', Form( "Tracking/NHitsSts%d", iDet), Form("NHitsSts%d_title", iDet),  100,0.,100.);
00091           fPosXSts_vs_time.at( iDet ) = MakeTH2('I', Form( "Tracking/XPosSts_vs_time%d", iDet), Form("XPosSts_vs_time%d_title", iDet),  3600,0,3600, 100,-50.,50.);
00092           fPosYSts_vs_time.at( iDet ) = MakeTH2('I', Form( "Tracking/YPosSts_vs_time%d", iDet), Form("YPosSts_vs_time%d_title", iDet),  3600,0,3600, 100,-50.,50.);
00093   }
00094 
00095 
00096 
00097   fPosXGem1 = MakeTH1('I', "Tracking/XPosGem1","XPosGem1_title",  100,-50.,50.);
00098   fPosYGem1 = MakeTH1('I', "Tracking/YPosGem1","YPosGem1_title",  100,-50.,50.);
00099   fPosXYGem1 = MakeTH2('I', "Tracking/XYPosGem1","XYPosGem1_title",  100,-50.,50.,  100,-50.,50.);
00100   fPosXGem2 = MakeTH1('I', "Tracking/XPosGem2","XPosGem2_title",  100,-50.,50.);
00101   fPosYGem2 = MakeTH1('I', "Tracking/YPosGem2","YPosGem2_title",  100,-50.,50.);
00102   fPosXYGem2 = MakeTH2('I', "Tracking/XYPosGem2","XYPosGem2_title",  100,-50.,50.,  100,-50.,50.);
00103 
00104   fResXSts1 = MakeTH1('I', "Tracking/XResSts1","XResSts1_title",  1000,-5.,5.);
00105   fResYSts1 = MakeTH1('I', "Tracking/YResSts1","YResSts1_title",  1000,-5.,5.);
00106   fResXSts2 = MakeTH1('I', "Tracking/XResSts2","XResSts2_title",  1000,-5.,5.);
00107   fResYSts2 = MakeTH1('I', "Tracking/YResSts2","YResSts2_title",  1000,-5.,5.);
00108   fResXSts3 = MakeTH1('I', "Tracking/XResSts3","XResSts3_title",  1000,-5.,5.);
00109   fResYSts3 = MakeTH1('I', "Tracking/YResSts3","YResSts3_title",  1000,-5.,5.);
00110 
00111   fResXGem1 = MakeTH1('I', "Tracking/XResGem1","XResGem1_title",  100,-5.,5.);
00112   fResYGem1 = MakeTH1('I', "Tracking/YResGem1","YResGem1_title",  100,-5.,5.);
00113   fResXGem2 = MakeTH1('I', "Tracking/XResGem2","XResGem2_title",  100,-5.,5.);
00114   fResYGem2 = MakeTH1('I', "Tracking/YResGem2","YResGem2_title",  100,-5.,5.);
00115 
00116  fPosXHod1_vs_time = MakeTH2('I', "Tracking/XPosHod1_vs_time","XPosHod1_vs_time_title",  3600,0,3600,100,-50.,50.);
00117  fPosYHod1_vs_time = MakeTH2('I', "Tracking/YPosHod1_vs_time","YPosHod1_vs_time_title",  3600,0,3600,100,-50.,50.);
00118  fPosXHod1_vs_time->GetXaxis()->SetTitle("time (sec) ");
00119  fPosXHod1_vs_time->GetYaxis()->SetTitle("H1_{hit} X (mm)");
00120  fPosYHod1_vs_time->GetXaxis()->SetTitle("time (sec) ");
00121  fPosYHod1_vs_time->GetYaxis()->SetTitle("H1_{hit} Y (mm)");
00122  fPosXHod2_vs_time = MakeTH2('I', "Tracking/XPosHod2_vs_time","XPosHod2_vs_time_title",  3600,0,3600,100,-50.,50.);
00123  fPosYHod2_vs_time = MakeTH2('I', "Tracking/YPosHod2_vs_time","YPosHod2_vs_time_title",  3600,0,3600,100,-50.,50.);
00124  fPosXHod2_vs_time->GetXaxis()->SetTitle("time (sec) ");
00125  fPosXHod2_vs_time->GetYaxis()->SetTitle("H2_{hit} X (mm)");
00126  fPosYHod2_vs_time->GetXaxis()->SetTitle("time (sec) ");
00127  fPosYHod2_vs_time->GetYaxis()->SetTitle("H2_{hit} Y (mm)");
00128 
00129 
00130 
00131 // JAM create full timestamp analysis histos if hodoscope provides full timestamp hits
00132 #ifdef HODOSCOPE_FULLTIMESTAMP
00133  fDeltaT_STS= MakeTH1('I',  Form( "Tracking/Times/DeltaT_Sts"), Form("Hit time difference within event for all sts"), 500000, -5000000., 20000., "#Delta t");
00134  fDeltaT_Hod= MakeTH1('I',  Form( "Tracking/Times/DeltaT_Hod"), Form("Hit time difference within event for all hodos"), 40000, -20000., 20000., "#Delta t");
00135  fDeltaT_GEM= MakeTH1('I',  Form( "Tracking/Times/DeltaT_GEM"), Form("Hit time difference within event for all gems"), 40000, -20000., 20000., "#Delta t");
00136  //fDeltaT_STS_HOD= MakeTH1('I',  Form( "Tracking/Times/DeltaT_STSvsHod"), Form("Hit time difference within event of all sts vs hod reference"), 5000, -5000., 5000., "#Delta t");
00137  fDeltaT_all= MakeTH1('I',  Form( "Tracking/Times/DeltaT_AllvsHod"), Form("Hit time difference within event of all components vs hod reference"), 40000, -20000., 20000., "#Delta t");
00138 #endif
00139 
00140 
00141 
00142 
00143    TString setupmacro = "set_TrackingPar.C";
00144    if (!gSystem->AccessPathName(setupmacro.Data())) {
00145       TGo4Log::Info("Executing tracking parameter setup script %s", setupmacro.Data());
00146       gROOT->ProcessLine(Form(".x %s", setupmacro.Data()));
00147   
00148    } else {
00149       TGo4Log::Info("NO tracking  setup script %s. Use previous values!", setupmacro.Data());
00150    }
00151 
00152 TString dir="Tracking";
00153 TString gname="HYGraph";
00154 TString fullname= dir+"/"+gname;
00155 cout <<"graph name is "<<fullname.Data() << endl;
00156  grHY=dynamic_cast<TGraphErrors*>(GetObject(fullname));
00157 
00158    if(grHY==0) {
00159       grHY=new TGraphErrors(2, ZHf, YHf, EZHf, EYHf);
00160       grHY->SetName(gname);
00161      AddObject(grHY, dir.Data());  
00162 
00163  }
00164  cout <<"assigned  object"<< grHY->GetName()<< endl;
00165 gname="HXGraph";
00166   fullname= dir+"/"+gname;
00167  grHX=dynamic_cast<TGraphErrors*>(GetObject(fullname));
00168    if(grHX==0) {
00169       grHX=new TGraphErrors(2, ZHf, XHf, EZHf, EXHf);
00170       grHX->SetName(gname);
00171       AddObject(grHX, dir.Data());  
00172 
00173  }
00174 gname="YGraph";
00175  fullname=  dir+"/"+gname;
00176   grY=dynamic_cast<TGraphErrors*>(GetObject(fullname));
00177   if(grY==0) {
00178         grY = new TGraphErrors(5, Zf, Yf, EZf, EYf);
00179         grY->SetName(gname); 
00180         grY->SetTitle("YZ projection");
00181         grY->GetXaxis()->SetTitle("Z (mm) ");
00182         grY->GetYaxis()->SetTitle("Y (mm) ");
00183         grY->SetMarkerStyle(21);
00184         grY->SetMarkerColor(4);
00185         grY->SetMarkerSize(0.8);
00186         grY->SetLineColor(4);
00187         grY->GetXaxis()->SetRangeUser(-10, 3100);
00188         AddObject(grY, dir.Data());  
00189  }
00190 gname="XGraph";
00191 fullname=  dir+"/"+gname;
00192  grX=dynamic_cast<TGraphErrors*>(GetObject(fullname));
00193    if(grX==0) {
00194         grX = new TGraphErrors(5, Zf, Xf, EZf, EXf);
00195         grX->SetName(gname);
00196         grX->SetTitle("XZ projection");
00197         grX->GetXaxis()->SetTitle("Z (mm) ");
00198         grX->GetYaxis()->SetTitle("X (mm) ");
00199         grX->SetMarkerStyle(20);
00200         grX->SetMarkerColor(4);
00201         grX->SetMarkerSize(0.8);
00202         grX->SetLineColor(4);
00203         grX->GetXaxis()->SetRangeUser(-10, 3100);
00204         AddObject(grX, dir.Data());  
00205  }
00206 gname="YaGraph";
00207  fullname=  dir+"/"+gname;
00208   graY=dynamic_cast<TGraphErrors*>(GetObject(fullname));
00209   if(graY==0) {
00210         graY = new TGraphErrors(7, Za, Ya, EZa, EYa);
00211         graY->SetName(gname); 
00212         graY->SetTitle("YZ projection");
00213         grY->GetXaxis()->SetTitle("Z (mm) ");
00214         graY->GetYaxis()->SetTitle("Y (mm) ");
00215         graY->SetMarkerStyle(21);
00216         graY->SetMarkerColor(4);
00217         graY->SetMarkerSize(0.8);
00218         graY->SetLineColor(4);
00219         graY->GetXaxis()->SetRangeUser(-10, 3100);
00220         AddObject(graY, dir.Data());  
00221  }
00222 gname="XaGraph";
00223 fullname=  dir+"/"+gname;
00224   graX=dynamic_cast<TGraphErrors*>(GetObject(fullname));
00225    if(graX==0) {
00226         graX = new TGraphErrors(7, Za, Xa, EZa, EXa);
00227         graX->SetName(gname);
00228         graX->SetTitle("XZ projection");
00229         graX->GetXaxis()->SetTitle("Z (mm) ");
00230         graX->GetYaxis()->SetTitle("X (mm) ");
00231         graX->SetMarkerStyle(20);
00232         graX->SetMarkerColor(4);
00233         graX->SetMarkerSize(0.8);
00234         graX->SetLineColor(4);
00235         graX->GetXaxis()->SetRangeUser(-10, 3100);
00236         AddObject(graX, dir.Data());  
00237  }
00238 
00239    printf("TTrackingProc Histograms created \n");
00240    fflush ( stdout);
00241 
00242 }
00243 
00244 
00245 
00246 TTrackingProc::~TTrackingProc()
00247 {
00248 }
00249 
00250 
00251 
00252 void TTrackingProc::InitEvent(TGo4EventElement* outevnt)
00253 {
00254    // first assign input event:
00255    // since input event object is never discarded within processor lifetime,
00256    // we just search for subevent by name once to speed up processing
00257  TCBMBeamtimeEvent* btevent=dynamic_cast<TCBMBeamtimeEvent*>(GetInputEvent());
00258       if(btevent)
00259       {
00260 //      fBeamMonitorEvent=dynamic_cast<TBeamMonitorEvent*>(btevent->GetSubEvent("BEAM"));
00261       fFirstHodoEvent=dynamic_cast<TFiberHodEvent*>(btevent->GetSubEvent("Hodo1"));
00262       fSecondHodoEvent=dynamic_cast<TFiberHodEvent*>(btevent->GetSubEvent("Hodo2"));
00263       fSTSEvent=dynamic_cast<TSTSEvent*>(btevent->GetSubEvent("STS"));
00264       fGEMEvent=dynamic_cast<TGEMEvent*>(btevent->GetSubEvent("GEM"));
00265       }
00266       else  {
00267          cout << "btevent == 0\n";
00268       }
00269    //if(fSTSEvent!=0)
00270  
00271  btevent=dynamic_cast<TCBMBeamtimeEvent*>(GetInputEvent("Detector"));
00272       if(btevent)
00273       {
00274         fRocInputEvent = dynamic_cast<TRocEvent*>(btevent->GetSubEvent("ROC"));
00275         fRocInputEvent->GetEventTime();
00276       }
00277       else  {
00278          cout << "detector event == 0\n";
00279       }
00280  
00281  
00282  
00283  // then assign output event
00284    // since output event object is never discarded within processor lifetime,
00285    // we just search for subevent by name once to speed up processing
00286    if(fOutputEvent==0)
00287    {
00288       TCBMBeamtimeEvent* btevent=dynamic_cast<TCBMBeamtimeEvent*>(outevnt); 
00289 // CbmKF* kalman = new CbmKF();
00290       if(btevent)
00291       {
00292          fOutputEvent=dynamic_cast<TTrackingEvent*>(btevent->GetSubEvent("TRACKING"));
00293       }
00294       else
00295       {
00296          fOutputEvent= dynamic_cast<TTrackingEvent*>(outevnt);
00297       }
00298       if(fOutputEvent==0) {
00299          GO4_STOP_ANALYSIS_MESSAGE("**** TTrackingProc: Fatal error: output event is not a TTrackingEvent!!! STOP GO4");
00300       }
00301    }
00302 }
00303 
00304 
00305  
00306 void TTrackingProc::FinalizeEvent()
00307 {
00308    if(fBeamMonitorEvent && fBeamMonitorEvent->IsValid())
00309    {
00310 
00311    }
00312    static int eventNum = 0;
00313    if( fFirstHodoEvent->NumHits() == 1 && fSecondHodoEvent->NumHits() == 1) {
00314 //initialisation of hits arrays
00315     eventNum++;
00316         for (int i2=0;i2<2;i2++)
00317         { 
00318           XHf[i2] = 0.0; 
00319           YHf[i2] = 0.0; 
00320         }
00321         for (int i5=0;i5<5;i5++)
00322         { 
00323           Xf[i5] = 0.0; 
00324           Yf[i5] = 0.0; 
00325         }
00326         for (int i7=0;i7<7;i7++)
00327         { 
00328           Xa[i7] = 0.0; 
00329           Ya[i7] = 0.0; 
00330         }
00331 
00332 #ifdef HODOSCOPE_FULLTIMESTAMP
00333          // initial time references for each system
00334          uint64_t HodoTime=0;
00335 //       uint64_t StsTime=0;
00336 //       uint64_t GemTime=0;
00337 #endif
00338 
00339 
00340 //      if(fFirstHodoEvent->NumHits()>0 && (fSecondHodoEvent->NumHits()+fGEMEvent->NumHits() > 0) ){
00341 //      if(fFirstHodoEvent->NumHits() > 0 &&  fSecondHodoEvent->NumHits() > 0 ){
00342 //              cout<<" H1 # of hits:  " << fFirstHodoEvent->NumHits()<<endl;
00343 //              cout<<" H2 # of hits:  " << fSecondHodoEvent->NumHits()<<endl;
00344 //              if(fGEMEvent) cout<<" GEM # of hits:  " << fGEMEvent->NumHits()<<endl;
00345 //              if(fSTSEvent) cout<<" Sts # of hits:  " << fSTSEvent->NumHits()<<endl;
00346 //      }
00347 //filling the control histos for the Hodo
00348         static double first_time_stamp = 0;
00349         if( eventNum == 1 )first_time_stamp = fRocInputEvent->GetEventTime();
00350         for( UInt_t iHit1 = 0; iHit1 < fFirstHodoEvent->NumHits(); iHit1++ ) {
00351             fPosXHod1->Fill( fFirstHodoEvent->Hit(iHit1).X );
00352             fPosYHod1->Fill( fFirstHodoEvent->Hit(iHit1).Y );
00353             fPosXYHod1->Fill( fFirstHodoEvent->Hit(iHit1).X , fFirstHodoEvent->Hit(iHit1).Y );
00354 //      cout<<"event # "<<eventNum<<" current time" << fRocInputEvent->GetEventTime()<<" first_time_stamp "<<first_time_stamp<<endl;
00355             //            fPosXHod1_vs_time->Fill( eventNum,fFirstHodoEvent->Hit(iHit1).X );
00356             fPosXHod1_vs_time->Fill( fRocInputEvent->GetEventTime()-first_time_stamp,fFirstHodoEvent->Hit(iHit1).X );
00357             fPosYHod1_vs_time->Fill( fRocInputEvent->GetEventTime()-first_time_stamp,fFirstHodoEvent->Hit(iHit1).Y );
00358             Xf[0] = fFirstHodoEvent->Hit(iHit1).X;
00359             Yf[0] = fFirstHodoEvent->Hit(iHit1).Y;
00360             Xa[0] = fFirstHodoEvent->Hit(iHit1).X;
00361             Ya[0] = fFirstHodoEvent->Hit(iHit1).Y;
00362 
00363 
00364 #ifdef HODOSCOPE_FULLTIMESTAMP
00365             uint64_t htime=fFirstHodoEvent->Hit(iHit1).fulltime;
00366             if(HodoTime==0)
00367                {
00368                   HodoTime=htime;
00369                }
00370             if(HodoTime)
00371                {
00372                   int64_t delta=htime-HodoTime;
00373                   fDeltaT_Hod->Fill(delta);
00374                   fDeltaT_all->Fill(delta);
00375                }
00376 #endif
00377 
00378         }
00379         for(UInt_t iHit5 = 0; iHit5 < fSecondHodoEvent->NumHits(); iHit5++ ) {
00380             fPosXHod2->Fill( fSecondHodoEvent->Hit(iHit5).X );
00381             fPosYHod2->Fill( fSecondHodoEvent->Hit(iHit5).Y );
00382             fPosXYHod2->Fill( fSecondHodoEvent->Hit(iHit5).X , fSecondHodoEvent->Hit(iHit5).Y ); 
00383             fPosXHod2_vs_time->Fill( fRocInputEvent->GetEventTime()-first_time_stamp,fSecondHodoEvent->Hit(iHit5).X );
00384             fPosYHod2_vs_time->Fill( fRocInputEvent->GetEventTime()-first_time_stamp,fSecondHodoEvent->Hit(iHit5).Y );
00385             Xf[4] = fSecondHodoEvent->Hit(iHit5).X;
00386             Yf[4] = fSecondHodoEvent->Hit(iHit5).Y;
00387             Xa[6] = fSecondHodoEvent->Hit(iHit5).X;
00388             Ya[6] = fSecondHodoEvent->Hit(iHit5).Y;
00389 
00390 #ifdef HODOSCOPE_FULLTIMESTAMP
00391             uint64_t htime=fSecondHodoEvent->Hit(iHit5).fulltime;
00392             if(HodoTime)
00393                      {
00394                         int64_t delta=htime-HodoTime;
00395                         fDeltaT_Hod->Fill(delta);
00396                         fDeltaT_all->Fill(delta);
00397                      }
00398 #endif
00399 
00400         }
00401 //    cout<<" Xa[0], Xa[6] " <<Xa[0]<<" "<< Xa[6]<<endl;
00402 //    cout<<" Ya[0], Ya[6] " <<Ya[0]<<" "<< Ya[6]<<endl;
00403             fCorrXHod->Fill( Xa[0] , Xa[6]); 
00404             fCorrYHod->Fill( Ya[0] , Ya[6]); 
00405 //          fCorrXHod->Fill( Xa[0] , Ya[6]); 
00406 //          fCorrYHod->Fill( Ya[0] , Xa[6]); 
00407         
00408 //filling the control histos for the Sts 
00409 
00410 
00411         for( Int_t iDet = 0; iDet < fSTSEvent->GetNPlanes(); iDet++ ) { //loop on stations
00412 //              fNHitsSts.at( iDet )->Fill( fSTSEvent->GetNPoints( iDet )); //loop on hits
00413                 for( Int_t iHit = 0; iHit < fSTSEvent->GetNPoints( iDet ); iHit++ ) {
00414                         fPosXSts.at( iDet )->Fill( -fSTSEvent->GetPoint( iDet, iHit ).GetX() );
00415                         if(iDet == 0) Xf[1] = fSTSEvent->GetPoint( iDet, iHit ).GetX();
00416                         if(iDet == 1) Xf[3] = fSTSEvent->GetPoint( iDet, iHit ).GetX();
00417                         if(iDet == 2) Xf[2] = fSTSEvent->GetPoint( iDet, iHit ).GetX();
00418                         fPosYSts.at( iDet )->Fill( fSTSEvent->GetPoint( iDet, iHit ).GetY() );
00419                         if(iDet == 0)Yf[1] = fSTSEvent->GetPoint( iDet, iHit ).GetY();
00420                         if(iDet == 1)Yf[3] = fSTSEvent->GetPoint( iDet, iHit ).GetY();
00421                         if(iDet == 2)Yf[2] = fSTSEvent->GetPoint( iDet, iHit ).GetY();
00422                         fPosXYSts.at( iDet )->Fill( fSTSEvent->GetPoint( iDet, iHit ).GetX(), fSTSEvent->GetPoint( iDet, iHit ).GetY() );
00423                         fPosXSts_vs_time.at( iDet )->Fill( fRocInputEvent->GetEventTime()-first_time_stamp, fSTSEvent->GetPoint( iDet, iHit ).GetX());
00424                         fPosYSts_vs_time.at( iDet )->Fill( fRocInputEvent->GetEventTime()-first_time_stamp, fSTSEvent->GetPoint( iDet, iHit ).GetY());
00425 
00426 #ifdef HODOSCOPE_FULLTIMESTAMP
00427          uint64_t stime=fSTSEvent->GetPoint( iDet, iHit ).GetTS();
00428          if(HodoTime && stime)
00429             {
00430                   static unsigned cnt=0;
00431                   int64_t delta=stime-HodoTime;
00432                   fDeltaT_STS->Fill(delta);
00433                   fDeltaT_all->Fill(delta);
00434 //                  if((delta > 10000 || delta < -10000) && (cnt++ % 100)==0)
00435 //                     {
00436 //                        cout <<"!!!!!! Too large sts to hodo time difference: dt="<<delta <<endl;
00437 //                        cout <<"!!!!!!         STStime="<<stime<<", hodotime="<<HodoTime <<endl;
00438 //
00439 //                     }
00440 
00441 
00442             }
00443 #endif
00444 
00445                 }
00446         }
00447 
00448 //filling the control histos for the GEMs
00449 //        for(UInt_t iHit6 = 0; iHit6 < fFirstHodoEvent->NumHits(); iHit6++ ) {
00450  int NumGem1Hits =0; 
00451  int NumGem2Hits =0; 
00452         for(UInt_t iHit6 = 0; iHit6 < fGEMEvent->NumHits(); iHit6++ ) {
00453             if(fGEMEvent->Hit(iHit6).StationIndex==0){
00454               fPosXGem1->Fill( -fGEMEvent->Hit(iHit6).X );
00455               fPosYGem1->Fill( -fGEMEvent->Hit(iHit6).Y );
00456               fPosXYGem1->Fill( -fGEMEvent->Hit(iHit6).X , -fGEMEvent->Hit(iHit6).Y );
00457               Xa[4] =  -fGEMEvent->Hit(iHit6).X; 
00458               Ya[4] =  -fGEMEvent->Hit(iHit6).Y;
00459 
00460 #ifdef HODOSCOPE_FULLTIMESTAMP
00461               uint64_t gtime=fGEMEvent->Hit(iHit6).Fulltime;
00462               if(HodoTime)
00463               {
00464                         int64_t delta=gtime-HodoTime;
00465                         fDeltaT_GEM->Fill(delta);
00466                         fDeltaT_all->Fill(delta);
00467               }
00468 #endif
00469 
00470               NumGem1Hits++;
00471             }
00472             if(fGEMEvent->Hit(iHit6).StationIndex==1){
00473               fPosXGem2->Fill( -fGEMEvent->Hit(iHit6).X );
00474               fPosYGem2->Fill( -fGEMEvent->Hit(iHit6).Y );
00475               fPosXYGem2->Fill( -fGEMEvent->Hit(iHit6).X , -fGEMEvent->Hit(iHit6).Y );
00476               Xa[5] =  -fGEMEvent->Hit(iHit6).X; 
00477               Ya[5] =  -fGEMEvent->Hit(iHit6).Y;
00478 
00479 #ifdef HODOSCOPE_FULLTIMESTAMP
00480               uint64_t gtime=fGEMEvent->Hit(iHit6).Fulltime;
00481         if(HodoTime)
00482              {
00483                   fDeltaT_GEM->Fill(gtime-HodoTime);
00484                   fDeltaT_all->Fill(gtime-HodoTime);
00485              }
00486 #endif
00487               NumGem2Hits++;
00488             }
00489         }
00490  //cout<<" event structure analysis" <<endl;
00491 //Sts event structure analysis 
00492         int NumSts1Hits = fSTSEvent->GetNPoints( 0); 
00493         int NumSts2Hits = fSTSEvent->GetNPoints( 2);  
00494         int NumSts3Hits = fSTSEvent->GetNPoints( 1); 
00495 //      NumSts1Hits =0; 
00496 //      NumSts2Hits = 0;  
00497 //      NumSts3Hits = 0; 
00498 if(NumSts1Hits==1 && NumSts3Hits==1)    fCorrXXSts.at(0)->Fill(Xf[1],Xf[3]);
00499 if(NumSts1Hits==1 && NumSts2Hits==1)    fCorrXXSts.at(1)->Fill(Xf[1],Xf[2]);
00500 if(NumSts2Hits==1 && NumSts3Hits==1)    fCorrXXSts.at(2)->Fill(Xf[2],Xf[3]);
00501 if(NumSts1Hits==1 && NumSts3Hits==1)    fCorrYYSts.at(0)->Fill(Yf[1],Yf[3]);
00502 if(NumSts1Hits==1 && NumSts2Hits==1)    fCorrYYSts.at(1)->Fill(Yf[1],Yf[2]);
00503 if(NumSts2Hits==1 && NumSts3Hits==1)    fCorrYYSts.at(2)->Fill(Yf[2],Yf[3]);
00504 
00505 
00506 if(NumSts1Hits==1)      fCorrXXH1Sts.at(0)->Fill(Xf[0],Xf[1]);
00507 if(NumSts2Hits==1)      fCorrXXH1Sts.at(1)->Fill(Xf[0],Xf[2]);
00508 if(NumSts3Hits==1)      fCorrXXH1Sts.at(2)->Fill(Xf[0],Xf[3]);
00509 if(NumSts1Hits==1)      fCorrYYH1Sts.at(0)->Fill(Yf[0],Yf[1]);
00510 if(NumSts2Hits==1)      fCorrYYH1Sts.at(1)->Fill(Yf[0],Yf[2]);
00511 if(NumSts3Hits==1)      fCorrYYH1Sts.at(2)->Fill(Yf[0],Yf[3]);
00512 
00513 if(Xf[0]>5. && Xf[0]<20. && Yf[0]>-2. && Yf[0]<5. &&Xf[4]>-2. && Xf[4]<17. && Yf[4]>7. && Yf[4]<22.){
00514         fNHitsSts.at(0)->Fill( NumSts1Hits); 
00515         fNHitsSts.at(1)->Fill( NumSts2Hits); 
00516         fNHitsSts.at(2)->Fill( NumSts3Hits); 
00517         }
00518 
00519 //H- Gem correlations
00520             if(NumGem1Hits == 1) {
00521             fCorrXHod1Gem1->Fill( Xa[0] , Xa[4]); 
00522             fCorrYHod1Gem1->Fill( Ya[0] , Ya[4]); 
00523             fCorrXHod2Gem1->Fill( Xa[6] , Xa[4]); 
00524             fCorrYHod2Gem1->Fill( Ya[6] , Ya[4]); 
00525            }
00526 
00527 //Sts- Gem correlations
00528             if(NumGem1Hits == 1&&NumSts1Hits==1) {
00529             fCorrXXGem1Sts.at(0)->Fill( Xf[1] , Xa[4]); 
00530             fCorrYYGem1Sts.at(0)->Fill( Yf[1] , Ya[4]); 
00531         }
00532             if(NumGem1Hits == 1&&NumSts2Hits==1) {
00533             fCorrXXGem1Sts.at(1)->Fill( Xf[2] , Xa[4]); 
00534             fCorrYYGem1Sts.at(1)->Fill( Yf[2] , Ya[4]);
00535         } 
00536             if(NumGem1Hits == 1&&NumSts3Hits==1) {
00537             fCorrXXGem1Sts.at(2)->Fill( Xf[3] , Xa[4]); 
00538             fCorrYYGem1Sts.at(2)->Fill( Yf[3] , Ya[4]);
00539         } 
00540 
00541 
00542 
00543 
00544 
00545 if(NumSts1Hits ==1&& NumSts2Hits ==1 && NumSts3Hits ==1) { //TRIGGER !!! all 5 points
00546 //if(NumSts1Hits ==1&& NumSts2Hits ==1 && NumSts3Hits ==1 && Xf[0]>5. && Xf[0]<20. && Yf[0]>-2. && Yf[0]<5. &&Xf[4]>-2. && Xf[4]<17. && Yf[4]>7. && Yf[4]<22. ) { //central (Sts) TRIGGER !!! all 5 points
00547 //if(NumSts1Hits ==1&& NumSts3Hits ==1) { //TRIGGER !!! STS1+STS3 +hodos
00548 //if(1) { //TRIGGER !!! hodos
00549 //if(NumSts1Hits == 1) { //TRIGGER !!! hodos +Sys1
00550                 if(Xf[1]<0.5&&Xf[1]>-0.5&&Yf[1]<0.5&&Yf[1]>-0.5){
00551 //              fPosXYHod1Sts->Fill(Xf[0],Yf[0]); H1
00552 //              fPosXYHod2Sts->Fill(Xf[4],Yf[4]); H2
00553                 fPosXYHod1Sts->Fill(Xf[2],Yf[2]);
00554                 fPosXYHod2Sts->Fill(Xf[3],Yf[3]);
00555                 }
00556 
00557 
00558 
00559         NPoint = 5; NHPoint = 2;
00560 //define Hodo alignment constants
00561         double ShiftX_Hod1 = 0.0;
00562         double ShiftY_Hod1 = 0.0;
00563         double ShiftX_Hod2 = 0.0;
00564         double ShiftY_Hod2 = 0.0;
00565 //define Sts alignment constants (mm)
00566         double ShiftX_Sts1 = 11.564; //run 155 (sig  = 2.6  mm)
00567         double ShiftY_Sts1 = 3.410;  //run 155 (sig  = 1.6 mm)
00568         double ShiftX_Sts2 =10.752;//run 155 (sig  = 2.9  mm)
00569         double ShiftY_Sts2 = 5.270; //run 155 (sig  = 1.9 mm)
00570         double ShiftX_Sts3 = 9.827;//run 155 (sig  = 3.4  mm)
00571         double ShiftY_Sts3 = 8.633; //run 155 (sig  = 2.3 mm)
00572 //define Gem alignment constants
00573         double ShiftX_Gem1 = -10.4; //mm
00574         double ShiftY_Gem1 = -11.5; //mm
00575         double ShiftX_Gem2 = 0.0;
00576         double ShiftY_Gem2 = 0.0;
00577 
00578         Zf[0] = 0.; Zf[1] = 270.; Zf[2] = 740.; Zf[3] = 1260.;  Zf[4] = 2770.; 
00579         Za[0] = 0.; Za[1] = 270.; Za[2] = 740.; Za[3] = 1260.;  Za[4] = 1630.; Za[5] = 1770.; Za[6] = 2770.; 
00580 // here should be Hodo and Sts hits reading part
00581         Xf[0] = Xf[0] + ShiftX_Hod1; 
00582         Xf[1] = Xf[1] + ShiftX_Sts1; 
00583         Xf[2] = Xf[2] + ShiftX_Sts2; 
00584         Xf[3] = Xf[3] + ShiftX_Sts3;  
00585         Xf[4] = Xf[4] + ShiftX_Hod2; 
00586         Yf[0] = Yf[0] + ShiftY_Hod1; 
00587         Yf[1] = Yf[1] + ShiftY_Sts1; 
00588         Yf[2] = Yf[2] + ShiftY_Sts2; 
00589         Yf[3] = Yf[3] + ShiftY_Sts3;  
00590         Yf[4] = Yf[4] + ShiftY_Hod2; 
00591 
00592         Xa[0] = Xa[0] + ShiftX_Hod1; 
00593         Xa[1] = Xf[1];
00594         Xa[2] = Xf[2]; 
00595         Xa[3] =Xf[3];  
00596         Xa[4] = Xa[4] + ShiftX_Gem1; 
00597         Xa[5] = Xa[5] + ShiftX_Gem2;
00598         Xa[6] = Xa[6] + ShiftX_Hod2; 
00599         
00600         Ya[0] = Ya[0] + ShiftY_Hod1; 
00601         Ya[1] = Yf[1]; 
00602         Ya[2] = Yf[2]; 
00603         Ya[3] = Yf[3];   
00604         Ya[4] = Ya[4] + ShiftY_Gem1; 
00605         Ya[5] = Ya[5] + ShiftY_Gem2; 
00606         Ya[6] = Ya[6] + ShiftY_Hod2; 
00607 //end of hits reading
00608 /*
00609         UInt_t HH1 = 0; UInt_t HH2 = 0;
00610         Xf[0] = fFirstHodoEvent->Hit(HH1).X;
00611         Xf[4] = fSecondHodoEvent->Hit(HH2).X;
00612         double XSts1 = (Zf[1] /Zf[4])*(Xf[4] - Xf[1]);
00613         double XSts2 = (Zf[2] /Zf[4])*(Xf[4] - Xf[2]);
00614         double XSts3 = (Zf[3] /Zf[4])*(Xf[4] - Xf[3]);
00615         double YSts1 = (Zf[1] /Zf[4])*(Yf[4] - Yf[1]);
00616         double YSts2 = (Zf[2] /Zf[4])*(Yf[4] - Yf[2]);
00617         double YSts3 = (Zf[3] /Zf[4])*(Yf[4] - Yf[3]);
00618         double distmin1 = 1000.;
00619         int smin1 =0;
00620         for (int s=0;s<NumSts1Hits;s++) 
00621         {
00622                 double dist = (XSts1 - fFirstStsEvent->Hit(s).X)*(XSts1 - fFirstStsEvent->Hit(s).X) + (YSts1 - fFirstStsEvent->Hit(s).Y)*(YSts1 - fFirstStsEvent->Hit(s).Y);
00623                 if(dist<distmin1){
00624                         distmin1 = dist;
00625                         smin1 = s;
00626                 }
00627         }
00628         double distmin2 = 1000.;
00629         int smin2 =0;
00630         for (int s=0;s<NumSts2Hits;s++) 
00631         {
00632                 double dist = (XSts2 - fSecondStsEvent->Hit(s).X)*(XSts2 - fSecondStsEvent->Hit(s).X) + (YSts2 - fSecondStsEvent->Hit(s).Y)*(YSts2 - fSecondStsEvent->Hit(s).Y);
00633                 if(dist<distmin2){
00634                         distmin2 = dist;
00635                         smin2 = s;
00636                 }
00637         }
00638         double distmin3 = 1000.;
00639         int smin3 =0;
00640         for (int s=0;s<NumSts3Hits;s++) 
00641         {
00642                 double dist = (XSts3 - fThirdStsEvent->Hit(s).X)*(XSts3 - fThirdStsEvent->Hit(s).X) + (YSts3 - fThirdStsEvent->Hit(s).Y)*(YSts3 - fThirdStsEvent->Hit(s).Y);
00643                 if(dist<distmin3){
00644                         distmin3 = dist;
00645                         smin3 = s;
00646                 }
00647         }
00648 */
00649 
00650 // define hits errors
00651         EZf[0] = 2.0; EZf[1] = 2.0; EZf[2] = 2.0; EZf[3] = 2.0;  EZf[4] = 2.0; 
00652         EXf[0] = 2.5; EXf[1] =0.014; EXf[2] = 0.014; EXf[3] = 0.014;  EXf[4] = 3.5; 
00653         EYf[0] = 2.0; EYf[1] =0.014; EYf[2] = 0.014; EYf[3] = 0.014;  EYf[4] = 2.8; 
00654 
00655         ZHf[0] =  Zf[0]; ZHf[1] = Zf[4]; 
00656         XHf[0] = Xf[0]; XHf[1] = Xf[4]; 
00657         YHf[0] = Yf[0]; YHf[1] = Yf[4]; 
00658         EZHf[0] = 2.0; EZHf[1] = 2.0; 
00659         EXHf[0] = 1.0;  EXHf[1] = 1.0; 
00660         EYHf[0] = 1.0; EYHf[1] = 1.0; 
00661 
00662         EZa[0] = 2.0; EZa[1] = 2.0; EZa[2] = 2.0; EZa[3] = 2.0;  EZa[4] = 2.0; EZa[5] = 2.0; EZa[6] = 2.0; 
00663         EXa[0] = 2.0; EXa[1] = 0.014; EXa[2] = 0.014; EXa[3] = 0.014;  EXa[4] = 1.5; EXa[5] = 1.5; EXa[6] = 2.0; 
00664         EYa[0] = 1.0; EYa[1] = 0.014; EYa[2] = 0.014; EYa[3] = 0.014;  EYa[4] = 0.5; EYa[5] = 0.5; EYa[6] = 1.0; 
00665 
00666         if(fFirstHodoEvent->NumHits()>0){
00667 //filling the Graphs
00668 for (int p=0;p<2;p++)
00669 { 
00670         //cout <<"p="<<p<<", grHY="<<(int )grHY<<", EZ="<<EZHf[p] <<", EY="<<EYHf[p]<< endl;
00671         grHY->SetPoint(p,ZHf[p], YHf[p]);
00672         grHY->SetPointError(p, EZHf[p], EYHf[p]);       
00673         grHX->SetPoint(p,ZHf[p], XHf[p]);
00674         grHX->SetPointError(p, EZHf[p], EXHf[p]);       
00675 }
00676 //cout<<" Hodo graphs are filled..."<<endl;
00677 
00678 for (int pp=0;pp<5;pp++)
00679 { 
00680         grY->SetPoint(pp,Zf[pp], Yf[pp]);
00681         grY->SetPointError(pp, EZf[pp], EYf[pp]);       
00682         grX->SetPoint(pp,Zf[pp], Xf[pp]);
00683         grX->SetPointError(pp, EZf[pp], EXf[pp]);       
00684 }
00685 //cout<<" Sts graphs are filled..."<<endl;
00686 
00687 for (int pp1=0;pp1<7;pp1++)
00688 { 
00689 //      cout<<" Za[pp1] "<<Za[pp1]<<"  Ya[pp1] "<<Ya[pp1]<<"  Xa[pp1] "<<Xa[pp1]<<endl;
00690         graY->SetPoint(pp1,Za[pp1], Ya[pp1]);
00691         graY->SetPointError(pp1, EZa[pp1], EYa[pp1]);   
00692         graX->SetPoint(pp1,Za[pp1], Xa[pp1]);
00693         graX->SetPointError(pp1, EZa[pp1], EXa[pp1]);   
00694 }
00695 
00696 //cout<<" All graphs are filled..."<<endl;
00697 if(NumGem2Hits <= 0) {
00698         graY->RemovePoint(5);
00699         graX->RemovePoint(5);
00700 }
00701 if(NumGem1Hits <= 0) {
00702         graY->RemovePoint(4);
00703         graX->RemovePoint(4);
00704 }
00705 
00706 if(NumSts3Hits <= 0) {
00707         grY->RemovePoint(3);
00708         grX->RemovePoint(3);
00709         graY->RemovePoint(3);
00710         graX->RemovePoint(3);
00711 }
00712 if(NumSts2Hits <= 0) {
00713         grY->RemovePoint(2);
00714         grX->RemovePoint(2);
00715         graY->RemovePoint(2);
00716         graX->RemovePoint(2);
00717 }
00718 if(NumSts1Hits <= 0) {
00719         grY->RemovePoint(1);
00720         grX->RemovePoint(1);
00721         graY->RemovePoint(1);
00722         graX->RemovePoint(1);
00723 }
00724 
00725 //fPosXYHod1Sts->Fill( Xf[0] , Yf[0] );
00726 //fPosXYHod2Sts->Fill( Xf[4] , Yf[4] );
00727 
00728 //fun = new TF1("fun", "[0]+[1]*x");
00729 //fun->SetLineColor(1);
00730 //grHY->Fit(fun,"Q");
00731 //grHY->Draw("line");
00732 //cout<<"Hodo Y Fit parameters "<<fun->GetParameter(0)<< " "<<fun->GetParameter(1)<<endl;
00733 //fun2 = new TF1("fun2", "[0]+[1]*x");
00734 //fun2->SetLineColor(2);
00735 //grHX->Draw("line");
00736 //cout<<"Hodo X Fit parameters "<<fun2->GetParameter(0)<< " "<<fun2->GetParameter(1)<<endl;
00737 fun3 = new TF1("fun3", "[0]+[1]*x"); 
00738 //fun3->SetLineColor(3);
00739 grY->Fit(fun3,"Q");
00740 //cout<<"All Y Fit parameters "<<fun3->GetParameter(0)<< " "<<fun3->GetParameter(1)<<endl;
00741 //cout<<"All Y Chi2/Ndf = "<<fun3->GetChisquare()<< "/ "<<fun3->GetNDF()<<"  = " <<fun3->GetChisquare()/fun3->GetNDF()<<endl;
00742 
00743 fun4 = new TF1("fun4", "[0]+[1]*x");
00744 //fun4->SetLineColor(4);
00745 grX->Fit(fun4,"Q");
00746 //cout<<"All X Fit parameters "<<fun4->GetParameter(0)<< " "<<fun4->GetParameter(1)<<endl;
00747 //cout<<"All X Chi2/Ndf = "<<fun4->GetChisquare()<< "/ "<<fun4->GetNDF()<<"  = " <<fun4->GetChisquare()/fun4->GetNDF()<<endl;
00748         //filling the Hodo's Residual histo:
00749 //            fResXHod1->Fill( XHf[0] - fun2->GetParameter(0) - fun2->GetParameter(1)*ZHf[0]);
00750 //            fResYHod1->Fill( YHf[0] - fun->GetParameter(0) - fun->GetParameter(1)*ZHf[0]);
00751 //            fResXHod2->Fill( XHf[1] - fun2->GetParameter(0) - fun2->GetParameter(1)*ZHf[1]);
00752 //            fResYHod2->Fill( YHf[1] - fun->GetParameter(0) - fun->GetParameter(1)*ZHf[1]);
00753 /*
00754             double ay = YHf[0];
00755             double by = (YHf[1]-YHf[0])/ZHf[1];
00756             double ax = XHf[0];
00757             double bx = (XHf[1]-XHf[0])/ZHf[1];
00758 */
00759 
00760 
00761 //            fResXHod1->Fill( XHf[0] - ax);
00762 //            fResYHod1->Fill( YHf[0] - ay);
00763 //            fResXHod2->Fill( XHf[1] - bx - ax*ZHf[1]);
00764 //            fResYHod2->Fill( YHf[1] - by - ay*ZHf[1]);
00765             fResXHod1->Fill( Xf[0] - fun4->GetParameter(0) - fun4->GetParameter(1)*Zf[0]);
00766             fResYHod1->Fill( Yf[0] - fun3->GetParameter(0) - fun3->GetParameter(1)*Zf[0]);
00767             fResXHod2->Fill( Xf[4] - fun4->GetParameter(0) - fun4->GetParameter(1)*Zf[4]);
00768             fResYHod2->Fill( Yf[4] - fun3->GetParameter(0) - fun3->GetParameter(1)*Zf[4]);
00769 //              fPosXYHod1Sts->Fill( fun4->GetParameter(0) + fun4->GetParameter(1)*Zf[0], fun3->GetParameter(0) + fun3->GetParameter(1)*Zf[0]);
00770 //              fPosXYHod2Sts->Fill( fun4->GetParameter(0) + fun4->GetParameter(1)*Zf[4], fun3->GetParameter(0) + fun3->GetParameter(1)*Zf[4]);
00771         //filling the Sts's Residual histo:
00772 if(NumSts1Hits > 0) {
00773             fResXSts1->Fill( Xf[1] - fun4->GetParameter(0) - fun4->GetParameter(1)*Zf[1]);
00774             fResYSts1->Fill( Yf[1] - fun3->GetParameter(0) - fun3->GetParameter(1)*Zf[1]);
00775         }
00776 if(NumSts2Hits > 0) {
00777             fResXSts2->Fill( Xf[2] - fun4->GetParameter(0) - fun4->GetParameter(1)*Zf[2]);
00778             fResYSts2->Fill( Yf[2] - fun3->GetParameter(0) - fun3->GetParameter(1)*Zf[2]);
00779         }
00780 if(NumSts3Hits > 0) {
00781             fResXSts3->Fill( Xf[3] - fun4->GetParameter(0) - fun4->GetParameter(1)*Zf[3]);
00782             fResYSts3->Fill( Yf[3] - fun3->GetParameter(0) - fun3->GetParameter(1)*Zf[3]);
00783         }
00784 fun5 = new TF1("fun5", "[0]+[1]*x"); 
00785 //fun3->SetLineColor(3);
00786 graY->Fit(fun5,"Q");
00787 //cout<<"All Ya Fit parameters "<<fun5->GetParameter(0)<< " "<<fun5->GetParameter(1)<<endl;
00788 //cout<<"All Ya Chi2/Ndf = "<<fun5->GetChisquare()<< "/ "<<fun5->GetNDF()<<"  = " <<fun5->GetChisquare()/fun5->GetNDF()<<endl;
00789 
00790 fun6 = new TF1("fun6", "[0]+[1]*x"); 
00791 graX->Fit(fun6,"Q");
00792 //cout<<"All Xa Fit parameters "<<fun6->GetParameter(0)<< " "<<fun6->GetParameter(1)<<endl;
00793 //cout<<"All Xa Chi2/Ndf = "<<fun6->GetChisquare()<< "/ "<<fun6->GetNDF()<<"  = " <<fun6->GetChisquare()/fun6->GetNDF()<<endl;
00794 
00795 if(NumGem1Hits > 0) {
00796             fResXGem1->Fill( Xa[4] - fun6->GetParameter(0) - fun6->GetParameter(1)*Za[4]);
00797             fResYGem1->Fill( Ya[4] - fun5->GetParameter(0) - fun5->GetParameter(1)*Za[4]);
00798         }
00799 if(NumGem2Hits > 0) {
00800             fResXGem2->Fill( Xa[5] - fun6->GetParameter(0) - fun6->GetParameter(1)*Za[5]);
00801             fResYGem2->Fill( Ya[5] - fun5->GetParameter(0) - fun5->GetParameter(1)*Za[5]);
00802         }
00803         
00804         
00805 //delete fun;
00806 //delete fun2;
00807 delete fun3;
00808 delete fun4;
00809 delete fun5;
00810 delete fun6;
00811                 }
00812         }
00813    }
00814    else {
00815 //      cout << " fFirstHodoEvent == 0 \n";
00816    }
00817 
00818 }