• Main Page
  • Related Pages
  • Namespaces
  • Data Structures
  • Files
  • File List
  • Globals

beamtime/tof-tdctest/go4/PLASTICS/TPlasticsProc.cxx (r4864/r4745)

Go to the documentation of this file.
00001 #include "TPlasticsProc.h"
00002 
00003 #include "TProfile.h"
00004 #include "TF1.h"
00005 #include "TMath.h"
00006 
00007 TPlasticsProc::TPlasticsProc(const char* name) : TCBMBeamtimeProc(name),
00008 fCrateInputEvent(0),
00009 fTriglogInputEvent(0),
00010 fVftxInputEvent(0),
00011 fGet4v1InputEvent(0),
00012 fOutputEvent(0)
00013 {
00014    TGo4Log::Info("**** TPlasticsProc: Create instance %s" , name );
00015 
00016    TString sName = name;
00017    UInt_t uIndexPlastic =((TString)sName(sName.Length()-2,2)).Atoi();
00018    fPar = (TPlasticsParam*) MakeParameter( Form("PlasticsPar_%02d", uIndexPlastic), "TPlasticsParam");
00019 
00020    Bool_t bAtLeastOneCaenPlastic = kFALSE;
00021    Bool_t bAtLeastOneVftxPlastic = kFALSE;
00022    Bool_t bAtLeastOneGet4Plastic = kFALSE;
00023    for( UInt_t uPlasticIndex = 0; uPlasticIndex < fPar->uNbPlastics; uPlasticIndex++)
00024       switch( fPar->uPlasticTdcType[uPlasticIndex] )
00025       {
00026          case 0:
00027             bAtLeastOneCaenPlastic = kTRUE;
00028             break;
00029          case 1:
00030             bAtLeastOneVftxPlastic = kTRUE;
00031             break;
00032          case 2:
00033             bAtLeastOneCaenPlastic = kTRUE;
00034             bAtLeastOneVftxPlastic = kTRUE;
00035             break;
00036          case 3:
00037             bAtLeastOneGet4Plastic = kTRUE;
00038             break;
00039       } // switch( fPar->uPlasticTdcType[uPlasticIndex] )
00040 
00041    if( kTRUE == bAtLeastOneVftxPlastic )
00042    {
00043       fVftxPar = (TVftxParam*) GetParameter("VftxPar");
00044       if( 0 == fVftxPar )
00045          fVftxPar = (TVftxParam*) MakeParameter("VftxPar", "TVftxParam");
00046    } // if( kTRUE == bAtLeastOneVftxPlastic )
00047 
00048    if( kTRUE == bAtLeastOneGet4Plastic)
00049    {
00050       // Try to recover the Get4 v1.0 Parameter object
00051       fGet4v1Par = (TGet4v1Param*) GetParameter("Get4v1Par");
00052       if( 0 == fGet4v1Par )
00053          fGet4v1Par = (TGet4v1Param*) MakeParameter("Get4v1Par", "TGet4v1Param");
00054    } // if( kTRUE == bAtLeastOneGet4Plastic)
00055 
00056    fParAnalysis = (TGsiAug12Param*) GetParameter("GsiAug12Par");
00057 
00058    //~ cout<<name<<endl;
00059 
00060    // To be put as preprocessor somewhere ? (T1290Data.h?)
00061    dCaenBinSize     = 25000.0/1024.0;
00062 
00063    if( 0 < fPar->uNbPlastics)
00064    {
00065       TString sFolder = Form("Plastics%02d", uIndexPlastic);
00066       fBeamProfilePlasticsTime = MakeTH2('D',
00067             Form("Plastics/%s/BeamProfilePlasticsTime_%02d", sFolder.Data(), uIndexPlastic),
00068             "Beam profile for Plastics" ,
00069             2400,-30000,30000,  // <= To be put as param later
00070             fPar->uNbPlastics, 0, fPar->uNbPlastics,
00071             "Left right dt [ps]", "Plastic index []" );
00072 
00073       Bool_t bOneCaenPlastic = kFALSE;
00074       Bool_t bOneGet4Plastic = kFALSE;
00075       for( UInt_t uPlasticIndex = 0; uPlasticIndex < fPar->uNbPlastics; uPlasticIndex++)
00076          if( 0 == fPar->uPlasticTdcType[uPlasticIndex] )
00077             bOneCaenPlastic = kTRUE;
00078             else if( 3 == fPar->uPlasticTdcType[uPlasticIndex] )
00079                bOneGet4Plastic = kTRUE;
00080       if( kTRUE == bOneCaenPlastic)
00081       {
00082          fReference2ProfilePlastics = MakeTH2('D',
00083                Form("Plastics/%s/DiamondProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00084                "Time profile against Diamond for Plastics" ,
00085                6000,-300000,300000,  // <= To be put as param later
00086                fPar->uNbPlastics, 0, fPar->uNbPlastics,
00087                "t(Diamond) - ( t(left)+ t(right) )/2 [ps]", "Plastic index []" );
00088 
00089          fReference1ProfilePlastics = MakeTH2('D',
00090                Form("Plastics/%s/ReferenceProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00091                "Time profile against Reference for Plastics" ,
00092                6000,-300000,300000,  // <= To be put as param later
00093                fPar->uNbPlastics, 0, fPar->uNbPlastics,
00094                "t(Reference) - ( t(left)+ t(right) )/2 [ps]", "Plastic index []" );
00095       } // if( kTRUE == bOneCaenPlastic)
00096          else if( kTRUE == bOneGet4Plastic )
00097          {
00098             fReference1ProfilePlastics = MakeTH2('D',
00099                   Form("Plastics/%s/Reference1ProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00100                   "Time profile against 1st Reference signal for Plastics" ,
00101                   6000,-300000,300000,  // <= To be put as param later
00102                   fPar->uNbPlastics, 0, fPar->uNbPlastics,
00103                   "t(Reference 1) - Mean Plastic Time [ps]", "Plastic index []" );
00104             fReference2ProfilePlastics = MakeTH2('D',
00105                   Form("Plastics/%s/Reference2ProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00106                   "Time profile against 2nd Reference signal (Diamond) for Plastics" ,
00107                   6000,-300000,300000,  // <= To be put as param later
00108                   fPar->uNbPlastics, 0, fPar->uNbPlastics,
00109                   "t(Reference 2)  - Mean Plastic Time [ps]", "Plastic index []" );
00110             fMeanRefProfilePlastics = MakeTH2('D',
00111                   Form("Plastics/%s/MeanRefProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00112                   "Profile against mean of 1st and 2nd reference signal for Plastics" ,
00113                   6000,-300000,300000,  // <= To be put as param later
00114                   fPar->uNbPlastics, 0, fPar->uNbPlastics,
00115                   "Distance Mean Time to Mean Reference [ps]", "Plastic index []" ) ;
00116          } // else if( kTRUE == bOneGet4Plastic )
00117          else
00118          {
00119             if( -1 < fVftxPar->iMainReferenceTdc ) // Only VFTX data!
00120             {
00121                if( -1 < fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc] )
00122                {
00123                   fReference1ProfilePlastics = MakeTH2('D',
00124                         Form("Plastics/%s/Reference1ProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00125                         "Time profile against 1st Reference signal for Plastics" ,
00126                         6000,-300000,300000,  // <= To be put as param later
00127                         fPar->uNbPlastics, 0, fPar->uNbPlastics,
00128                         "t(Reference 1) - Mean Plastic Time [ps]", "Plastic index []" );
00129                   if( -1 < fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc] )
00130                      fMeanRefProfilePlastics = MakeTH2('D',
00131                            Form("Plastics/%s/MeanRefProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00132                            "Profile against mean of 1st and 2nd reference signal for Plastics" ,
00133                            6000,-300000,300000,  // <= To be put as param later
00134                            fPar->uNbPlastics, 0, fPar->uNbPlastics,
00135                            "Distance Mean Time to Mean Reference [ps]", "Plastic index []" ) ;
00136                } // if( -1 < fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc] )
00137                if( -1 < fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc] )
00138                   fReference2ProfilePlastics = MakeTH2('D',
00139                         Form("Plastics/%s/Reference2ProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00140                         "Time profile against 2nd Reference signal (Diamond) for Plastics" ,
00141                         6000,-300000,300000,  // <= To be put as param later
00142                         fPar->uNbPlastics, 0, fPar->uNbPlastics,
00143                         "t(Reference 2)  - Mean Plastic Time [ps]", "Plastic index []" );
00144             } // else if( -1 < fVftxPar->iMainReferenceTdc ) // Contains VFTX data!
00145          }
00146 
00147       // Maybe put the size of thez Tot histo in time as an option in parameter file
00148       fTotLeftPlastics = MakeTH2('D',
00149             Form("Plastics/%s/TotLeftProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00150             "Left Tot profile for Plastics" ,
00151             2000, -50, 50,
00152             fPar->uNbPlastics, 0. , (Double_t)fPar->uNbPlastics,
00153             "Left Tot [ns]", "Plastic index []" ) ; ;
00154 
00155       // Maybe put the size of thez Tot histo in time as an option in parameter file
00156       fTotRightPlastics = MakeTH2('D',
00157             Form("Plastics/%s/TotRightProfilePlastics_%02d", sFolder.Data(), uIndexPlastic),
00158             "Right Tot profile for Plastics" ,
00159             2000, -50, 50,
00160             fPar->uNbPlastics, 0. , (Double_t)fPar->uNbPlastics,
00161             "Right Tot [ns]", "Plastics index []" ) ;
00162 
00163       // Multiplicity plots
00164       fMultiplicityPlastics = MakeTH2('I',
00165             Form("Plastics/%s/MultiplicityPlastics_%02d", sFolder.Data(), uIndexPlastic),
00166             "Plastics Multiplicity " ,
00167             fPar->uNbPlastics, 0. , (Double_t)fPar->uNbPlastics,
00168             TVftxBoardData::MaxMult, 0, TVftxBoardData::MaxMult,
00169             "Plastics index []",
00170             "Multiplicity (Number of Hits/event) []" ) ;
00171    } // if( 0 < fPar->uNbPlastics)
00172 
00173    TGo4Log::Info("**** TPlasticsProc: Instance %s created", name );
00174 }
00175 
00176 
00177 TPlasticsProc::~TPlasticsProc()
00178 {
00179    cout << "**** TPlasticsProc: Delete instance " << endl;
00180 }
00181 
00182 void TPlasticsProc::InitEvent(TGo4EventElement* outevnt)
00183 {
00184    // first assign input event:
00185    // since input event object is never discarded within processor lifetime,
00186    // we just search for subevent by name once to speed up processing
00187    Bool_t bAtLeastOneCaenPlastic = kFALSE;
00188    Bool_t bAtLeastOneVftxPlastic = kFALSE;
00189    Bool_t bAtLeastOneGet4Plastic = kFALSE;
00190    for( UInt_t uPlasticIndex = 0; uPlasticIndex < fPar->uNbPlastics; uPlasticIndex++)
00191       switch( fPar->uPlasticTdcType[uPlasticIndex] )
00192       {
00193          case 0:
00194             bAtLeastOneCaenPlastic = kTRUE;
00195             break;
00196          case 1:
00197             bAtLeastOneVftxPlastic = kTRUE;
00198             break;
00199          case 2:
00200             bAtLeastOneCaenPlastic = kTRUE;
00201             bAtLeastOneVftxPlastic = kTRUE;
00202             break;
00203          case 3:
00204             bAtLeastOneGet4Plastic = kTRUE;
00205             break;
00206       } // switch( fPar->uPlasticTdcType[uPlasticIndex] )
00207 
00208 
00209    if(fCrateInputEvent==0 && kTRUE == bAtLeastOneCaenPlastic)
00210    {
00211      TCBMBeamtimeEvent* btevent=dynamic_cast<TCBMBeamtimeEvent*>(GetOutputEvent("Unpack"));
00212           if(btevent)
00213           {
00214                  fCrateInputEvent=dynamic_cast<TMbsCrateEvent*>(btevent->GetSubEvent("MBSCRATE"));
00215           }
00216           else
00217           {
00218                  fCrateInputEvent=dynamic_cast<TMbsCrateEvent*>(GetInputEvent());
00219           }
00220           if(fCrateInputEvent==0) {
00221        GO4_STOP_ANALYSIS_MESSAGE("**** TPlasticsProc: Fatal error: Unpack output event has no TMbsCrateEvent!!! STOP GO4");
00222           }
00223    } // if(fCrateInputEvent==0 && kTRUE == bAtLeastOneCaenPlastic)
00224    if(fTriglogInputEvent==0 && kTRUE == fParAnalysis->bWithTriglog)
00225    {
00226      TCBMBeamtimeEvent* btevent=dynamic_cast<TCBMBeamtimeEvent*>(GetOutputEvent("Unpack"));
00227           if(btevent)
00228           {
00229                  fTriglogInputEvent=dynamic_cast<TTriglogEvent*>(btevent->GetSubEvent("TRIGLOG"));
00230           }
00231 
00232           if(fTriglogInputEvent==0) {
00233        GO4_STOP_ANALYSIS_MESSAGE("**** TRpcProc: Fatal error: Unpack output event is/has not a TTriglogEvent!!! STOP GO4");
00234           }
00235    } // if(fTriglogInputEvent==0)
00236    if(fVftxInputEvent==0 && kTRUE == bAtLeastOneVftxPlastic)
00237    {
00238           TCBMBeamtimeEvent* btevent=dynamic_cast<TCBMBeamtimeEvent*>(GetInputEvent());
00239           if(btevent)
00240           {
00241                  fVftxInputEvent=dynamic_cast<TVftxEvent*>(btevent->GetSubEvent("VFTX"));
00242           }
00243 
00244           if(fVftxInputEvent==0) {
00245                  GO4_STOP_ANALYSIS_MESSAGE("**** TPlasticsProc: Fatal error: input event has no TVftxEvent!!! STOP GO4");
00246           }
00247    } // if(fVftxInputEvent==0 && kTRUE == bAtLeastOneVftxPlastic)
00248    if(fGet4v1InputEvent==0 && kTRUE == bAtLeastOneGet4Plastic)
00249    {
00250      TCBMBeamtimeEvent* btevent=dynamic_cast<TCBMBeamtimeEvent*>(GetOutputEvent("Unpack"));
00251      if(btevent)
00252      {
00253         fGet4v1InputEvent=dynamic_cast<TGet4v1Event*>(btevent->GetSubEvent("ROCGET4V10"));
00254      }
00255 
00256      if(fGet4v1InputEvent==0) {
00257         GO4_STOP_ANALYSIS_MESSAGE("**** TPlasticsProc: Fatal error: Unpack output event has no TGet4v1Event!!! STOP GO4");
00258      }
00259    } // if(fGet4v1InputEvent==0 && kTRUE == bAtLeastOneGet4Plastic)
00260 
00261    // then assign output event
00262    // since output event object is never discarded within processor lifetime,
00263    // we just search for subevent by name once to speed up processing
00264    if(fOutputEvent==0)
00265    {
00266       TCBMBeamtimeEvent* btevent=dynamic_cast<TCBMBeamtimeEvent*>(outevnt);
00267       if(btevent)
00268       {
00269          TString sName = GetName();
00270          if ( sName.Contains("PLASTICS_") )
00271          {
00272             fOutputEvent=dynamic_cast<TPlasticsEvent*>(btevent->GetSubEvent(
00273                   (TString)(sName(sName.Index("PLASTICS_"), sName.Length() ) ) ) );
00274          }
00275          else
00276             fOutputEvent=dynamic_cast<TPlasticsEvent*>(btevent->GetSubEvent( "PLASTICS" ) );
00277       }
00278       else
00279       {
00280          fOutputEvent= dynamic_cast<TPlasticsEvent*>(outevnt);
00281       }
00282       if(fOutputEvent==0) {
00283          GO4_STOP_ANALYSIS_MESSAGE("**** TPlasticsProc: Fatal error: output event is not a TPlasticsEvent!!! STOP GO4");
00284       }
00285       else {
00286          //     BuildEvent(dynamic_cast<TGo4MbsSubEvent*>(btevent->GetSubEvent("MBSCRATE")));
00287       }
00288    }
00289 
00290 }
00291 
00292 void TPlasticsProc::FinalizeEvent()
00293 {
00294    hitCurrent.Clear();
00295 
00296    if( -1 < fParAnalysis->iTriggerRejection && kTRUE == fParAnalysis->bWithTriglog)
00297       if( 1 == ((fTriglogInputEvent->fVulomTriggerPattern>>fParAnalysis->iTriggerRejection)&0x1 )  )
00298       {
00299           // The trigger Rejection is active and we have rejected trigger
00300          return;
00301       }
00302 
00303    // No plastic is defined => skip processing
00304    if( 0 == fPar->uNbPlastics )
00305    {
00306       fOutputEvent->SetValid(kTRUE);
00307       return;
00308    }
00309 
00310    for( UInt_t uPlasticIndex = 0; uPlasticIndex < fPar->uNbPlastics; uPlasticIndex++)
00311    {
00312       if( 3 == fPar->uPlasticTdcType[uPlasticIndex] )
00313       {
00314          if( 0 == fGet4v1InputEvent->fEvents.size() )
00315          {
00316             fMultiplicityPlastics->Fill(uPlasticIndex, 0 );
00317             if( kFALSE == fGet4v1Par->bFreeStreaming || kTRUE == fGet4v1Par->bRawDataMode )
00318             {
00319                if( 0 == (fOutputEvent->fEvents).size() )
00320                      (fOutputEvent->fEvents).resize( 1 );
00321                ((fOutputEvent->fEvents)[0].fHits[uPlasticIndex]).clear();
00322             } // if( kFALSE == fGet4v1Par->bFreeStreaming || kTRUE == fGet4v1Par->bRawDataMode )
00323             continue;
00324          } // if( 0 == fGet4v1InputEvent->fEvents.size() )
00325          if( (fOutputEvent->fEvents).size() != fGet4v1InputEvent->fEvents.size() )
00326             (fOutputEvent->fEvents).resize( fGet4v1InputEvent->fEvents.size() );
00327 
00328          for( UInt_t uGet4EventIndex = 0; uGet4EventIndex < fGet4v1InputEvent->fEvents.size(); uGet4EventIndex++)
00329          {
00330             if( 2 == fPar->uNbSides[uPlasticIndex] )
00331                ProcessGet4v10PlasticDouble( uGet4EventIndex, uPlasticIndex );
00332                else ProcessGet4v10PlasticSingle( uGet4EventIndex, uPlasticIndex );
00333 
00334             fMultiplicityPlastics->Fill(uPlasticIndex,
00335                   ((fOutputEvent->fEvents[uGet4EventIndex]).fHits[uPlasticIndex]).size() );
00336          } // for( UInt_t uGet4EventIndex = 0; uGet4EventIndex < fGet4v1InputEvent->fEvents.size(); uGet4EventIndex++)
00337       } // if( 3 == fPar->fPar->uPlasticTdcType[uPlasticIndex] ) => GET4 v1.0 !!!
00338       else if(  fPar->uPlasticTdcType[uPlasticIndex] < 3 )
00339       {
00340          (fOutputEvent->fEvents).resize(1);
00341          if( 1 == fPar->uPlasticTdcType[uPlasticIndex] )
00342          {
00343             // The trigger selection is active and we have another trigger => skip processing
00344             if( kTRUE == fParAnalysis->bWithTriglog )
00345                if( !( ( -1 == fPar->iTriggerSelection) ||
00346                    ( 1 == ((fTriglogInputEvent->fVulomTriggerPattern>>fPar->iTriggerSelection)&0x1 )  ) ) )
00347                   continue;
00348 
00349             // Check that the VFTX data are valid after calibration
00350             if(fVftxInputEvent && fVftxInputEvent->IsValid())
00351             {
00352                if( 2 == fPar->uNbSides[uPlasticIndex] )
00353                   ProcessVftxPlasticDouble( uPlasticIndex );
00354                   else ProcessVftxPlasticSingle( uPlasticIndex );
00355             } // if(fVftxInputEvent && fVftxInputEvent->IsValid())
00356          } // else if( 1 == fPar->uPlasticTdcType[uPlasticIndex] ) => VFTX  !!!
00357          else
00358          {
00359             // The trigger selection is active and we have another trigger => skip processing
00360             if( !( ( -1 == fPar->iTriggerSelection) ||
00361                 ( 1 == ((fTriglogInputEvent->fVulomTriggerPattern>>fPar->iTriggerSelection)&0x1 )  ) ) )
00362                continue;
00363             // Check that the MBS unpack event is valid
00364             if(fCrateInputEvent && fCrateInputEvent->IsValid())
00365             {
00366                if( 2 == fPar->uNbSides[uPlasticIndex] )
00367                   ProcessCaenPlasticDouble( uPlasticIndex );
00368                   else ProcessCaenPlasticSingle( uPlasticIndex );
00369             } // if(fCrateInputEvent && fCrateInputEvent->IsValid())
00370 
00371          } // else of if( 1 == fPar->uPlasticTdcType[uPlasticIndex] ) => CAEN  !!!
00372          fMultiplicityPlastics->Fill(uPlasticIndex,
00373                ((fOutputEvent->fEvents[0]).fHits[uPlasticIndex]).size() );
00374       } // else if( fPar->uPlasticTdcType[uPlasticIndex] < 3 ) => NOT Get4 v1.0
00375    } // for( UInt_t uPlasticIndex = 0; uPlasticIndex < fPar->uNbPlastics; uPlasticIndex++)
00376 
00377    fOutputEvent->SetValid(kTRUE);
00378 }
00379 
00380 /**************************************************************************************************/
00381 void TPlasticsProc::ProcessCaenPlasticSingle( UInt_t uPlasticIndex, Double_t dCaenOtherOffset )
00382 {
00383    // Multiple TDC data check
00384    if( 1 < (fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).hit_lead[
00385                                                      fPar->uChannelPlasticLeftTime[uPlasticIndex] ] ||
00386         1 < (fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTot[uPlasticIndex]  ]).hit_trail[
00387                                                              fPar->uChannelPlasticLeftTot[uPlasticIndex] ])
00388    {
00389       (fOutputEvent->fEvents[0]).fbMultiEdgesPresent = kTRUE;
00390       hitCurrent.fbMultiEdge = kTRUE;
00391    } // If at least one edge/side has more than 1 TDC data
00392 
00393    // Test if both time are there
00394    for(Int_t iHitIndex = 0;
00395          iHitIndex<(fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).hit_lead[
00396                                                      fPar->uChannelPlasticLeftTime[uPlasticIndex] ] ;
00397          iHitIndex++)
00398    {
00399       // Test if both time are there
00400       if( -1 < (fCrateInputEvent->fMtdc[
00401                   fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00402                            fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex] )
00403       {
00404          if( -1 < (fCrateInputEvent->fMtdc[
00405                            fPar->uTdcPlasticLeftTot[uPlasticIndex]  ]).trail_multi[
00406                               fPar->uChannelPlasticLeftTot[uPlasticIndex] ][iHitIndex] )
00407          {
00408             fTotLeftPlastics->Fill(
00409                ( (fCrateInputEvent->fMtdc[
00410                      fPar->uTdcPlasticLeftTot[uPlasticIndex] ]).trail_multi[
00411                         fPar->uChannelPlasticLeftTot[uPlasticIndex]][iHitIndex]
00412                  -(fCrateInputEvent->fMtdc[
00413                      fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00414                         fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex]
00415                )*dCaenBinSize*fPar->dToTGainListLeft[uPlasticIndex]/1000.0
00416                - fPar->dOffsetListLeft[uPlasticIndex]/1000.0
00417                - fPar->dTotOffsetListLeft[uPlasticIndex],
00418                (Double_t)uPlasticIndex  );
00419 
00420             if( -1 < fPar->iDiamondTdcVFTX && -1 < fPar->iDiamondChannelVFTX )
00421                if( -1 < (fCrateInputEvent->fMtdc[
00422                             fPar->iDiamondTdcVFTX ]).lead_multi[
00423                             fPar->iDiamondChannelVFTX][0] )
00424             {
00425               fReference2ProfilePlastics->Fill(
00426                ( ( (fCrateInputEvent->fMtdc[fPar->iDiamondTdcVFTX ]).lead_multi[
00427                       fPar->iDiamondChannelVFTX][0])*dCaenBinSize -
00428                    ( (  (fCrateInputEvent->fMtdc[
00429                          fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00430                             fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex]
00431                      )*dCaenBinSize
00432                      + fPar->dOffsetListLeft[uPlasticIndex]
00433                    )
00434                  ),
00435                 (Double_t)uPlasticIndex );
00436             } // if Diamond info defined and Diamond time there
00437             if( -1 < fPar->iReferenceTdc && -1 < fPar->iReferenceChannel )
00438                if( -1 < (fCrateInputEvent->fMtdc[
00439                             fPar->iReferenceTdc ]).lead_multi[
00440                             fPar->iReferenceChannel][0])
00441             {
00442                fReference1ProfilePlastics->Fill(
00443                 ( ( (fCrateInputEvent->fMtdc[fPar->iReferenceTdc ]).lead_multi[
00444                        fPar->iReferenceChannel][0])*dCaenBinSize -
00445                    ( ( (fCrateInputEvent->fMtdc[
00446                          fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00447                             fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex]
00448                      )*dCaenBinSize
00449                      + fPar->dOffsetListLeft[uPlasticIndex]
00450                    )
00451                  ),
00452                 (Double_t)uPlasticIndex );
00453             } // if Reference info defined and Reference time there
00454 
00455             // Store time including offset!
00456             hitCurrent.dTimeLeft  = dCaenBinSize*( (fCrateInputEvent->fMtdc[
00457                                       fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00458                                            fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex] )
00459                                     + fPar->dOffsetListLeft[uPlasticIndex]
00460                                     - dCaenOtherOffset;
00461             // Store tot including offset and gain!
00462             // Missing the offset for the coarse counter overflow!!!!
00463             hitCurrent.dTotLeft   =
00464                ( (fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTot[uPlasticIndex] ]).trail_multi[
00465                                        fPar->uChannelPlasticLeftTot[uPlasticIndex]][iHitIndex]
00466                 -(fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00467                                        fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex]
00468                )*dCaenBinSize*fPar->dToTGainListLeft[uPlasticIndex]
00469                - fPar->dOffsetListLeft[uPlasticIndex]
00470                - fPar->dTotOffsetListLeft[uPlasticIndex]*1000.0;
00471 
00472             ((fOutputEvent->fEvents[0]).fHits[uPlasticIndex]).push_back( hitCurrent );
00473 
00474             hitCurrent.Clear();
00475          } // if valid Tot on left side
00476       } // if( valid on both ends!)
00477    } // for TdcDataIndex in both end
00478    return;
00479 }
00480 void TPlasticsProc::ProcessCaenPlasticDouble( UInt_t uPlasticIndex, Double_t dCaenOtherOffset )
00481 {
00482    // Multiple TDC data check
00483    if( 1 < (fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).hit_lead[
00484                                                      fPar->uChannelPlasticLeftTime[uPlasticIndex] ] ||
00485        1 < (fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTot[uPlasticIndex]  ]).hit_trail[
00486                                                      fPar->uChannelPlasticLeftTot[uPlasticIndex] ] ||
00487        1 < (fCrateInputEvent->fMtdc[ fPar->uTdcPlasticRightTime[uPlasticIndex]  ]).hit_lead[
00488                                                      fPar->uChannelPlasticRightTime[uPlasticIndex] ] ||
00489        1 < (fCrateInputEvent->fMtdc[ fPar->uTdcPlasticRightTot[uPlasticIndex]  ]).hit_trail[
00490                                                      fPar->uChannelPlasticRightTot[uPlasticIndex] ])
00491    {
00492       (fOutputEvent->fEvents[0]).fbMultiEdgesPresent = kTRUE;
00493       hitCurrent.fbMultiEdge = kTRUE;
00494    } // If at least one edge/side has more than 1 TDC data
00495 
00496    // Test if both time are there
00497    for(Int_t iHitIndex = 0;
00498          iHitIndex<(fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).hit_lead[
00499                                                      fPar->uChannelPlasticLeftTime[uPlasticIndex] ] &&
00500           iHitIndex<(fCrateInputEvent->fMtdc[ fPar->uTdcPlasticRightTime[uPlasticIndex]  ]).hit_lead[
00501                                                      fPar->uChannelPlasticRightTime[uPlasticIndex] ];
00502          iHitIndex++)
00503    {
00504       // Test if both time are there
00505       if( -1 < (fCrateInputEvent->fMtdc[
00506                   fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00507                            fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex] &&
00508           -1 < (fCrateInputEvent->fMtdc[
00509                         fPar->uTdcPlasticRightTime[uPlasticIndex] ]).lead_multi[
00510                            fPar->uChannelPlasticRightTime[uPlasticIndex]][iHitIndex] )
00511       {
00512          fBeamProfilePlasticsTime->Fill(
00513                                ( ( (fCrateInputEvent->fMtdc[
00514                                        fPar->uTdcPlasticRightTime[uPlasticIndex] ]).lead_multi[
00515                                           fPar->uChannelPlasticRightTime[uPlasticIndex]][iHitIndex]
00516                                   -(fCrateInputEvent->fMtdc[
00517                                        fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00518                                           fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex]
00519                                 )*dCaenBinSize
00520                                 + fPar->dOffsetListRight[uPlasticIndex]
00521                                 - fPar->dOffsetListLeft[uPlasticIndex]
00522                                      ),
00523                ((Double_t)uPlasticIndex) +0.00001 );
00524 
00525          if( -1 < (fCrateInputEvent->fMtdc[
00526                            fPar->uTdcPlasticRightTot[uPlasticIndex]  ]).trail_multi[
00527                               fPar->uChannelPlasticRightTot[uPlasticIndex] ][iHitIndex] )
00528          {
00529             fTotRightPlastics->Fill(
00530                ( (fCrateInputEvent->fMtdc[
00531                      fPar->uTdcPlasticRightTot[uPlasticIndex] ]).trail_multi[
00532                         fPar->uChannelPlasticRightTot[uPlasticIndex]][iHitIndex]
00533                  -(fCrateInputEvent->fMtdc[
00534                      fPar->uTdcPlasticRightTime[uPlasticIndex]  ]).lead_multi[
00535                         fPar->uChannelPlasticRightTime[uPlasticIndex] ][iHitIndex]
00536                )*dCaenBinSize*fPar->dToTGainListRight[uPlasticIndex]/1000.0
00537                - fPar->dOffsetListRight[uPlasticIndex]/1000.0
00538                - fPar->dTotOffsetListRight[uPlasticIndex],
00539                (Double_t)uPlasticIndex  );
00540          } // if valid Tot on right side
00541 
00542          if( -1 < (fCrateInputEvent->fMtdc[
00543                            fPar->uTdcPlasticLeftTot[uPlasticIndex]  ]).trail_multi[
00544                               fPar->uChannelPlasticLeftTot[uPlasticIndex] ][iHitIndex] )
00545          {
00546             fTotLeftPlastics->Fill(
00547                ( (fCrateInputEvent->fMtdc[
00548                      fPar->uTdcPlasticLeftTot[uPlasticIndex] ]).trail_multi[
00549                         fPar->uChannelPlasticLeftTot[uPlasticIndex]][iHitIndex]
00550                  -(fCrateInputEvent->fMtdc[
00551                      fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00552                         fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex]
00553                )*dCaenBinSize*fPar->dToTGainListLeft[uPlasticIndex]/1000.0
00554                - fPar->dOffsetListLeft[uPlasticIndex]/1000.0
00555                - fPar->dTotOffsetListLeft[uPlasticIndex],
00556                (Double_t)uPlasticIndex  );
00557 
00558             if( -1 < (fCrateInputEvent->fMtdc[
00559                               fPar->uTdcPlasticRightTot[uPlasticIndex]  ]).trail_multi[
00560                                  fPar->uChannelPlasticRightTot[uPlasticIndex] ][iHitIndex] )
00561             {
00562                if( -1 < fPar->iDiamondTdcVFTX && -1 < fPar->iDiamondChannelVFTX )
00563                   if( -1 < (fCrateInputEvent->fMtdc[
00564                                fPar->iDiamondTdcVFTX ]).lead_multi[
00565                                fPar->iDiamondChannelVFTX][0] )
00566                {
00567                  fReference2ProfilePlastics->Fill(
00568                   ( ( ( (fCrateInputEvent->fMtdc[fPar->iDiamondTdcVFTX ]).lead_multi[
00569                          fPar->iDiamondChannelVFTX][0])*dCaenBinSize -
00570                       ( (  (fCrateInputEvent->fMtdc[
00571                             fPar->uTdcPlasticRightTime[uPlasticIndex] ]).lead_multi[
00572                                fPar->uChannelPlasticRightTime[uPlasticIndex]][iHitIndex]
00573                          + (fCrateInputEvent->fMtdc[
00574                             fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00575                                fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex]
00576                         )*dCaenBinSize
00577                         + fPar->dOffsetListRight[uPlasticIndex]
00578                         + fPar->dOffsetListLeft[uPlasticIndex]
00579                       )/2.
00580                           ) ),
00581                    (Double_t)uPlasticIndex );
00582                } // if Diamond info defined and Diamond time there
00583                if( -1 < fPar->iReferenceTdc && -1 < fPar->iReferenceChannel )
00584                   if( -1 < (fCrateInputEvent->fMtdc[
00585                                fPar->iReferenceTdc ]).lead_multi[
00586                                fPar->iReferenceChannel][0])
00587                {
00588                   fReference1ProfilePlastics->Fill(
00589                   ( ( ( (fCrateInputEvent->fMtdc[fPar->iReferenceTdc ]).lead_multi[
00590                           fPar->iReferenceChannel][0])*dCaenBinSize -
00591                       ( (  (fCrateInputEvent->fMtdc[
00592                             fPar->uTdcPlasticRightTime[uPlasticIndex] ]).lead_multi[
00593                                fPar->uChannelPlasticRightTime[uPlasticIndex]][iHitIndex]
00594                          + (fCrateInputEvent->fMtdc[
00595                             fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00596                                fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex]
00597                         )*dCaenBinSize
00598                         + fPar->dOffsetListRight[uPlasticIndex]
00599                         + fPar->dOffsetListLeft[uPlasticIndex]
00600                       )/2.
00601                           ) ),
00602                    (Double_t)uPlasticIndex );
00603                } // if Reference info defined and Reference time there
00604 
00605                // Store time including offset!
00606                hitCurrent.dTimeLeft  = dCaenBinSize*( (fCrateInputEvent->fMtdc[
00607                                          fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00608                                               fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex] )
00609                                        + fPar->dOffsetListLeft[uPlasticIndex]
00610                                        - dCaenOtherOffset;
00611                // Store tot including offset and gain!
00612                // Missing the offset for the coarse counter overflow!!!!
00613                hitCurrent.dTotLeft   =
00614                   ( (fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTot[uPlasticIndex] ]).trail_multi[
00615                                           fPar->uChannelPlasticLeftTot[uPlasticIndex]][iHitIndex]
00616                    -(fCrateInputEvent->fMtdc[ fPar->uTdcPlasticLeftTime[uPlasticIndex]  ]).lead_multi[
00617                                           fPar->uChannelPlasticLeftTime[uPlasticIndex] ][iHitIndex]
00618                   )*dCaenBinSize*fPar->dToTGainListLeft[uPlasticIndex]
00619                   - fPar->dOffsetListLeft[uPlasticIndex]
00620                   - fPar->dTotOffsetListLeft[uPlasticIndex]*1000.0;
00621                // Store time including offset!
00622                hitCurrent.dTimeRight = dCaenBinSize*( (fCrateInputEvent->fMtdc[
00623                                            fPar->uTdcPlasticRightTime[uPlasticIndex] ]).lead_multi[
00624                                               fPar->uChannelPlasticRightTime[uPlasticIndex]][iHitIndex] )
00625                                        + fPar->dOffsetListRight[uPlasticIndex]
00626                                        - dCaenOtherOffset;
00627                // Store tot including offset and gain!
00628                // Missing the offset for the coarse counter overflow!!!!
00629                hitCurrent.dTotRight  =
00630                   ( (fCrateInputEvent->fMtdc[ fPar->uTdcPlasticRightTot[uPlasticIndex] ]).trail_multi[
00631                                           fPar->uChannelPlasticRightTot[uPlasticIndex]][iHitIndex]
00632                    -(fCrateInputEvent->fMtdc[ fPar->uTdcPlasticRightTime[uPlasticIndex]  ]).lead_multi[
00633                                           fPar->uChannelPlasticRightTime[uPlasticIndex] ][iHitIndex]
00634                   )*dCaenBinSize*fPar->dToTGainListRight[uPlasticIndex]
00635                   - fPar->dOffsetListRight[uPlasticIndex]
00636                   - fPar->dTotOffsetListRight[uPlasticIndex]*1000.0;
00637 
00638                ((fOutputEvent->fEvents[0]).fHits[uPlasticIndex]).push_back( hitCurrent );
00639 
00640                hitCurrent.Clear();
00641             }// if valid trailing edge also on right side => valid hit!
00642          } // if valid Tot on left side
00643       } // if( valid on both ends!)
00644    } // for TdcDataIndex in both end
00645    return;
00646 }
00647 /**************************************************************************************************/
00648 void TPlasticsProc::ProcessVftxPlasticSingle( UInt_t uPlasticIndex, Double_t dVftxOtherOffset )
00649 {
00650    // Multiple TDC data check
00651    if( 1 < (fVftxInputEvent->fVftxBoards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).iMultiplicity[
00652                                                      fPar->uChannelPlasticLeftTime[uPlasticIndex] ]  ||
00653        1 < (fVftxInputEvent->fVftxBoards[ fPar->uTdcPlasticLeftTot[uPlasticIndex] ] ).iMultiplicity[
00654                                                      fPar->uChannelPlasticLeftTot[uPlasticIndex] ]  )
00655    {
00656       (fOutputEvent->fEvents[0]).fbMultiEdgesPresent = kTRUE;
00657       hitCurrent.fbMultiEdge = kTRUE;
00658    } // If at least one edge/side has more than 1 TDC data
00659 
00660    // Too close Hit merging attempt
00661    UInt_t uTotTdcHitToUseLeft  = 0;
00662 
00663    // Test if both time are there
00664    for(Int_t iHitIndex = 0;
00665          iHitIndex < (fVftxInputEvent->fVftxBoards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).iMultiplicity[
00666                    fPar->uChannelPlasticLeftTime[uPlasticIndex] ];
00667          iHitIndex++)
00668    {
00669       // Too close Hit merging attempt
00670       uTotTdcHitToUseLeft  = iHitIndex;
00671 
00672       if( kTRUE == fVftxInputEvent->IsHitThere(
00673                   fPar->uTdcPlasticLeftTime[uPlasticIndex],
00674                   fPar->uChannelPlasticLeftTime[uPlasticIndex],
00675                   iHitIndex ) )
00676       {
00677          Double_t dAutoOffset = CLOCK_TIME*(
00678                fVftxPar->iAutomaticTdcOffset[ fPar->uTdcPlasticLeftTime[ uPlasticIndex] ] );
00679 
00680          // Too close Hit merging attempt
00681          if( 0 < fPar->dMinimalTimeBetweenHits &&
00682              iHitIndex+1 < (fVftxInputEvent->fVftxBoards[
00683                             fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).iMultiplicity[
00684                             fPar->uChannelPlasticLeftTime[uPlasticIndex] ]  &&
00685              iHitIndex+1 < (fVftxInputEvent->fVftxBoards[
00686                             fPar->uTdcPlasticLeftTot[uPlasticIndex] ] ).iMultiplicity[
00687                             fPar->uChannelPlasticLeftTot[uPlasticIndex] ] &&
00688              kTRUE == fVftxInputEvent->IsHitThere(
00689                          fPar->uTdcPlasticLeftTot[uPlasticIndex],
00690                          fPar->uChannelPlasticLeftTot[uPlasticIndex],
00691                          iHitIndex ) &&
00692              kTRUE == fVftxInputEvent->IsHitThere(
00693                          fPar->uTdcPlasticLeftTot[uPlasticIndex],
00694                          fPar->uChannelPlasticLeftTot[uPlasticIndex],
00695                          iHitIndex+1 ) &&
00696              kTRUE == fVftxInputEvent->IsHitThere(
00697                         fPar->uTdcPlasticLeftTime[uPlasticIndex],
00698                         fPar->uChannelPlasticLeftTime[uPlasticIndex],
00699                         iHitIndex+ 1 ) )
00700          {
00701             Double_t dDistanceBtwnFirstSecond =
00702                   fVftxInputEvent->GetCalibratedTime(
00703                        fPar->uTdcPlasticLeftTime[uPlasticIndex],
00704                        fPar->uChannelPlasticLeftTime[uPlasticIndex],
00705                        iHitIndex+ 1, fVftxPar->uUseCoarseCorrectedTime )
00706                 - fVftxInputEvent->GetCalibratedTime(
00707                        fPar->uTdcPlasticLeftTot[uPlasticIndex],
00708                        fPar->uChannelPlasticLeftTot[uPlasticIndex],
00709                        iHitIndex, fVftxPar->uUseCoarseCorrectedTime );
00710             // For now take merge all hits under limit, even those
00711             // where second hit starts befor first hit finish
00712             if( 0 < dDistanceBtwnFirstSecond &&
00713                 dDistanceBtwnFirstSecond < fPar->dMinimalTimeBetweenHits )
00714                uTotTdcHitToUseLeft = iHitIndex+ 1;
00715          } // if min hit dist defined and at least 2 full hits
00716 
00717          // Check Left tot is there
00718          if( kTRUE == fVftxInputEvent->IsHitThere(
00719                         fPar->uTdcPlasticLeftTot[uPlasticIndex],
00720                         fPar->uChannelPlasticLeftTot[uPlasticIndex],
00721                         iHitIndex ) )
00722          {
00723             fTotLeftPlastics->Fill(
00724                ( fVftxInputEvent->GetCalibratedTime(
00725                      fPar->uTdcPlasticLeftTot[uPlasticIndex],
00726                      fPar->uChannelPlasticLeftTot[uPlasticIndex],
00727                      uTotTdcHitToUseLeft, fVftxPar->uUseCoarseCorrectedTime )
00728                  - fVftxInputEvent->GetCalibratedTime(
00729                        fPar->uTdcPlasticLeftTime[uPlasticIndex],
00730                        fPar->uChannelPlasticLeftTime[uPlasticIndex],
00731                        iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
00732                 )*fPar->dToTGainListLeft[uPlasticIndex]/1000.0
00733 //                - fPar->dOffsetListLeft[uPlasticIndex]/1000.0
00734                 - fPar->dTotOffsetListLeft[uPlasticIndex],
00735                (Double_t)uPlasticIndex );
00736 
00737             // Check first if the main reference TDC is defined in VFTX parameter
00738             if( -1 < fVftxPar->iMainReferenceTdc )
00739             {
00740                // Then Check if the 1st reference channel is defined for this TDC
00741                if( -1 < fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc] )
00742                   if( kTRUE == fVftxInputEvent->IsHitThere(
00743                                  fVftxPar->iMainReferenceTdc,
00744                                  fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc],
00745                                  0 ) )
00746                   {
00747                      fReference1ProfilePlastics->Fill(
00748                      ( fVftxInputEvent->GetCalibratedTime(
00749                            fVftxPar->iMainReferenceTdc,
00750                            fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc],
00751                            0, fVftxPar->uUseCoarseCorrectedTime )
00752                        - ( fVftxInputEvent->GetCalibratedTime(
00753                                fPar->uTdcPlasticLeftTime[uPlasticIndex],
00754                                fPar->uChannelPlasticLeftTime[uPlasticIndex],
00755                                iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
00756                            + fPar->dOffsetListLeft[uPlasticIndex]
00757                            - ( 1 == fVftxPar->uTdcOffsetEnable? dAutoOffset : 0.0)
00758                          ) ),
00759                       (Double_t)uPlasticIndex );
00760                      // Then Check if the 2nd reference channel is defined for this TDC
00761                      if( -1 < fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc] )
00762                         if( kTRUE == fVftxInputEvent->IsHitThere(
00763                                        fVftxPar->iMainReferenceTdc,
00764                                        fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc],
00765                                        0 ) )
00766                         {
00767                            Double_t dRefMeanSum = 0.;
00768                            dRefMeanSum += fVftxInputEvent->GetCalibratedTime(
00769                                              fVftxPar->iMainReferenceTdc,
00770                                              fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc],
00771                                              0, fVftxPar->uUseCoarseCorrectedTime );
00772                            dRefMeanSum += fVftxInputEvent->GetCalibratedTime(
00773                                              fVftxPar->iMainReferenceTdc,
00774                                              fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc],
00775                                              0, fVftxPar->uUseCoarseCorrectedTime );
00776                            fMeanRefProfilePlastics->Fill(
00777                             ( dRefMeanSum/2. -
00778                                 ( fVftxInputEvent->GetCalibratedTime(
00779                                       fPar->uTdcPlasticLeftTime[uPlasticIndex],
00780                                       fPar->uChannelPlasticLeftTime[uPlasticIndex],
00781                                       iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
00782                                   + fPar->dOffsetListLeft[uPlasticIndex]
00783                                   - ( 1 == fVftxPar->uTdcOffsetEnable? dAutoOffset : 0.0)
00784                                 ) ),
00785                              (Double_t)uPlasticIndex );
00786                         }
00787                   }
00788                // Then Check if the 2nd reference channel is defined for this TDC
00789                if( -1 < fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc] )
00790                   if( kTRUE == fVftxInputEvent->IsHitThere(
00791                                  fVftxPar->iMainReferenceTdc,
00792                                  fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc],
00793                                  0 ) )
00794                   {
00795                      fReference2ProfilePlastics->Fill(
00796                       ( fVftxInputEvent->GetCalibratedTime(
00797                             fVftxPar->iMainReferenceTdc,
00798                             fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc],
00799                             0, fVftxPar->uUseCoarseCorrectedTime )
00800                         - ( fVftxInputEvent->GetCalibratedTime(
00801                                 fPar->uTdcPlasticLeftTime[uPlasticIndex],
00802                                 fPar->uChannelPlasticLeftTime[uPlasticIndex],
00803                                 iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
00804                             + fPar->dOffsetListLeft[uPlasticIndex]
00805                             - ( 1 == fVftxPar->uTdcOffsetEnable? dAutoOffset : 0.0)
00806                           ) ),
00807                        (Double_t)uPlasticIndex );
00808                   }
00809             } // if( -1 < fVftxPar->iMainReferenceTdc )
00810 
00811             // Store time including offset!
00812             hitCurrent.dTimeLeft  = fVftxInputEvent->GetCalibratedTime(
00813                                           fPar->uTdcPlasticLeftTime[uPlasticIndex],
00814                                           fPar->uChannelPlasticLeftTime[uPlasticIndex],
00815                                           iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
00816                                     + fPar->dOffsetListLeft[uPlasticIndex]
00817                                     - dVftxOtherOffset
00818                                     - CLOCK_TIME*( 1 == fVftxPar->uTdcOffsetEnable?
00819                                           fVftxPar->iAutomaticTdcOffset[ fPar->uTdcPlasticLeftTime[ uPlasticIndex] ] : 0.0);
00820 
00821             // Store tot including offset and gain!
00822             // Missing the offset for the coarse counter overflow!!!!
00823             hitCurrent.dTotLeft   =
00824                ( fVftxInputEvent->GetCalibratedTime(
00825                      fPar->uTdcPlasticLeftTot[uPlasticIndex],
00826                      fPar->uChannelPlasticLeftTot[uPlasticIndex],
00827                      uTotTdcHitToUseLeft, fVftxPar->uUseCoarseCorrectedTime )
00828                 - fVftxInputEvent->GetCalibratedTime(
00829                       fPar->uTdcPlasticLeftTime[uPlasticIndex],
00830                       fPar->uChannelPlasticLeftTime[uPlasticIndex],
00831                       iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
00832                )*fPar->dToTGainListLeft[uPlasticIndex]
00833 //               - fPar->dOffsetListLeft[uPlasticIndex]
00834                - fPar->dTotOffsetListLeft[uPlasticIndex]*1000.0;
00835 
00836             ((fOutputEvent->fEvents[0]).fHits[uPlasticIndex]).push_back( hitCurrent );
00837 
00838             hitCurrent.Clear();
00839          } // if valid trailing edge on left side
00840       } // if( valid on both ends!)
00841 
00842       // Too close Hit merging attempt
00843       if( iHitIndex + 1 == (Int_t)uTotTdcHitToUseLeft )
00844          iHitIndex = uTotTdcHitToUseLeft;
00845    } // for TdcDataIndex in both end
00846 
00847    return;
00848 }
00849 void TPlasticsProc::ProcessVftxPlasticDouble( UInt_t uPlasticIndex, Double_t dVftxOtherOffset )
00850 {
00851    // Multiple TDC data check
00852    if( 1 < (fVftxInputEvent->fVftxBoards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).iMultiplicity[
00853                                                      fPar->uChannelPlasticLeftTime[uPlasticIndex] ]  ||
00854        1 < (fVftxInputEvent->fVftxBoards[ fPar->uTdcPlasticLeftTot[uPlasticIndex] ] ).iMultiplicity[
00855                                                      fPar->uChannelPlasticLeftTot[uPlasticIndex] ]   ||
00856        1 < (fVftxInputEvent->fVftxBoards[ fPar->uTdcPlasticRightTime[uPlasticIndex] ]  ).iMultiplicity[
00857                                                      fPar->uChannelPlasticRightTime[uPlasticIndex] ] ||
00858        1 < (fVftxInputEvent->fVftxBoards[ fPar->uTdcPlasticRightTot[uPlasticIndex] ] ).iMultiplicity[
00859                                                      fPar->uChannelPlasticRightTot[uPlasticIndex] ] )
00860    {
00861       (fOutputEvent->fEvents[0]).fbMultiEdgesPresent = kTRUE;
00862       hitCurrent.fbMultiEdge = kTRUE;
00863    } // If at least one edge/side has more than 1 TDC data
00864 
00865    // Too close Hit merging attempt
00866    UInt_t uTotTdcHitToUseLeft  = 0;
00867    UInt_t uTotTdcHitToUseRight = 0;
00868 
00869    // Test if both time are there
00870    for(Int_t iHitIndex = 0;
00871          iHitIndex< (fVftxInputEvent->fVftxBoards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).iMultiplicity[
00872                    fPar->uChannelPlasticLeftTime[uPlasticIndex] ] &&
00873          iHitIndex< (fVftxInputEvent->fVftxBoards[ fPar->uTdcPlasticRightTime[uPlasticIndex] ]  ).iMultiplicity[
00874                     fPar->uChannelPlasticRightTime[uPlasticIndex] ];
00875          iHitIndex++)
00876    {
00877       // Too close Hit merging attempt
00878       uTotTdcHitToUseLeft  = iHitIndex;
00879       uTotTdcHitToUseRight = iHitIndex;
00880 
00881       if( kTRUE == fVftxInputEvent->IsHitThere(
00882                      fPar->uTdcPlasticLeftTime[uPlasticIndex],
00883                      fPar->uChannelPlasticLeftTime[uPlasticIndex],
00884                      iHitIndex )  &&
00885           kTRUE == fVftxInputEvent->IsHitThere(
00886                      fPar->uTdcPlasticRightTime[uPlasticIndex],
00887                      fPar->uChannelPlasticRightTime[uPlasticIndex],
00888                      iHitIndex ) )
00889       {
00890          Double_t dAutoOffset = CLOCK_TIME*(
00891                fVftxPar->iAutomaticTdcOffset[ fPar->uTdcPlasticRightTime[uPlasticIndex] ]
00892               -fVftxPar->iAutomaticTdcOffset[ fPar->uTdcPlasticLeftTime[ uPlasticIndex] ] );
00893          Double_t dAutoOffsetSum = CLOCK_TIME*(
00894                fVftxPar->iAutomaticTdcOffset[ fPar->uTdcPlasticRightTime[uPlasticIndex] ]
00895               +fVftxPar->iAutomaticTdcOffset[ fPar->uTdcPlasticLeftTime[ uPlasticIndex] ] );
00896 
00897          fBeamProfilePlasticsTime->Fill(
00898                                ( fVftxInputEvent->GetCalibratedTime(
00899                                      fPar->uTdcPlasticRightTime[uPlasticIndex],
00900                                      fPar->uChannelPlasticRightTime[uPlasticIndex],
00901                                      iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
00902                                 - fVftxInputEvent->GetCalibratedTime(
00903                                       fPar->uTdcPlasticLeftTime[uPlasticIndex],
00904                                       fPar->uChannelPlasticLeftTime[uPlasticIndex],
00905                                       iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
00906                                 + fPar->dOffsetListRight[uPlasticIndex]
00907                                 - fPar->dOffsetListLeft[uPlasticIndex]
00908                                 - ( 1 == fVftxPar->uTdcOffsetEnable? dAutoOffset : 0.0)
00909                                      ),
00910                ((Double_t)uPlasticIndex) +0.00001 );
00911 
00912          // Too close Hit merging attempt
00913          if( 0 < fPar->dMinimalTimeBetweenHits &&
00914              iHitIndex+1 < (fVftxInputEvent->fVftxBoards[
00915                             fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).iMultiplicity[
00916                             fPar->uChannelPlasticLeftTime[uPlasticIndex] ]  &&
00917              iHitIndex+1 < (fVftxInputEvent->fVftxBoards[
00918                             fPar->uTdcPlasticLeftTot[uPlasticIndex] ] ).iMultiplicity[
00919                             fPar->uChannelPlasticLeftTot[uPlasticIndex] ] &&
00920              kTRUE == fVftxInputEvent->IsHitThere(
00921                          fPar->uTdcPlasticLeftTot[uPlasticIndex],
00922                          fPar->uChannelPlasticLeftTot[uPlasticIndex],
00923                          iHitIndex ) &&
00924              kTRUE == fVftxInputEvent->IsHitThere(
00925                          fPar->uTdcPlasticLeftTot[uPlasticIndex],
00926                          fPar->uChannelPlasticLeftTot[uPlasticIndex],
00927                          iHitIndex+1 ) &&
00928              kTRUE == fVftxInputEvent->IsHitThere(
00929                         fPar->uTdcPlasticLeftTime[uPlasticIndex],
00930                         fPar->uChannelPlasticLeftTime[uPlasticIndex],
00931                         iHitIndex+ 1 ) )
00932          {
00933             Double_t dDistanceBtwnFirstSecond =
00934                   fVftxInputEvent->GetCalibratedTime(
00935                        fPar->uTdcPlasticLeftTime[uPlasticIndex],
00936                        fPar->uChannelPlasticLeftTime[uPlasticIndex],
00937                        iHitIndex+ 1, fVftxPar->uUseCoarseCorrectedTime )
00938                 - fVftxInputEvent->GetCalibratedTime(
00939                        fPar->uTdcPlasticLeftTot[uPlasticIndex],
00940                        fPar->uChannelPlasticLeftTot[uPlasticIndex],
00941                        iHitIndex, fVftxPar->uUseCoarseCorrectedTime );
00942             // For now take merge all hits under limit, even those
00943             // where second hit starts befor first hit finish
00944             if( 0 < dDistanceBtwnFirstSecond &&
00945                 dDistanceBtwnFirstSecond < fPar->dMinimalTimeBetweenHits )
00946                uTotTdcHitToUseLeft = iHitIndex+ 1;
00947          } // if min hit dist defined and at least 2 full hits
00948          if( 0 < fPar->dMinimalTimeBetweenHits &&
00949              iHitIndex+1 < (fVftxInputEvent->fVftxBoards[
00950                             fPar->uTdcPlasticRightTime[uPlasticIndex] ] ).iMultiplicity[
00951                             fPar->uChannelPlasticRightTime[uPlasticIndex] ]  &&
00952              iHitIndex+1 < (fVftxInputEvent->fVftxBoards[
00953                             fPar->uTdcPlasticRightTot[uPlasticIndex] ] ).iMultiplicity[
00954                             fPar->uChannelPlasticRightTot[uPlasticIndex] ] &&
00955              kTRUE == fVftxInputEvent->IsHitThere(
00956                          fPar->uTdcPlasticRightTot[uPlasticIndex],
00957                          fPar->uChannelPlasticRightTot[uPlasticIndex],
00958                          iHitIndex ) &&
00959              kTRUE == fVftxInputEvent->IsHitThere(
00960                          fPar->uTdcPlasticRightTot[uPlasticIndex],
00961                          fPar->uChannelPlasticRightTot[uPlasticIndex],
00962                          iHitIndex+1 ) &&
00963              kTRUE == fVftxInputEvent->IsHitThere(
00964                         fPar->uTdcPlasticRightTime[uPlasticIndex],
00965                         fPar->uChannelPlasticRightTime[uPlasticIndex],
00966                         iHitIndex+ 1 ) )
00967          {
00968             Double_t dDistanceBtwnFirstSecond =
00969                   fVftxInputEvent->GetCalibratedTime(
00970                        fPar->uTdcPlasticRightTime[uPlasticIndex],
00971                        fPar->uChannelPlasticRightTime[uPlasticIndex],
00972                        iHitIndex+ 1, fVftxPar->uUseCoarseCorrectedTime )
00973                 - fVftxInputEvent->GetCalibratedTime(
00974                        fPar->uTdcPlasticRightTot[uPlasticIndex],
00975                        fPar->uChannelPlasticRightTot[uPlasticIndex],
00976                        iHitIndex, fVftxPar->uUseCoarseCorrectedTime );
00977             // For now take merge all hits under limit, even those
00978             // where second hit starts befor first hit finish
00979             if( 0 < dDistanceBtwnFirstSecond &&
00980                 dDistanceBtwnFirstSecond < fPar->dMinimalTimeBetweenHits )
00981                uTotTdcHitToUseRight = iHitIndex+ 1;
00982          } // if min hit dist defined and at least 2 full hits
00983 
00984          // Check Right tot is there
00985          if( kTRUE == fVftxInputEvent->IsHitThere(
00986                      fPar->uTdcPlasticRightTot[uPlasticIndex],
00987                      fPar->uChannelPlasticRightTot[uPlasticIndex],
00988                      uTotTdcHitToUseRight ) )
00989          {
00990             fTotRightPlastics->Fill(
00991                ( fVftxInputEvent->GetCalibratedTime(
00992                        fPar->uTdcPlasticRightTot[uPlasticIndex],
00993                        fPar->uChannelPlasticRightTot[uPlasticIndex],
00994                        uTotTdcHitToUseRight, fVftxPar->uUseCoarseCorrectedTime )
00995                  - fVftxInputEvent->GetCalibratedTime(
00996                        fPar->uTdcPlasticRightTime[uPlasticIndex],
00997                        fPar->uChannelPlasticRightTime[uPlasticIndex],
00998                        iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
00999                 )*fPar->dToTGainListRight[uPlasticIndex]/1000.0
01000 //                - fPar->dOffsetListRight[uPlasticIndex]/1000.0
01001                 - fPar->dTotOffsetListRight[uPlasticIndex],
01002                (Double_t)uPlasticIndex );
01003          } // if valid trailing edge on right side
01004          // Check Left tot is there
01005          if( kTRUE == fVftxInputEvent->IsHitThere(
01006                      fPar->uTdcPlasticLeftTot[uPlasticIndex],
01007                      fPar->uChannelPlasticLeftTot[uPlasticIndex],
01008                      uTotTdcHitToUseLeft ) )
01009          {
01010             fTotLeftPlastics->Fill(
01011                ( fVftxInputEvent->GetCalibratedTime(
01012                      fPar->uTdcPlasticLeftTot[uPlasticIndex],
01013                      fPar->uChannelPlasticLeftTot[uPlasticIndex],
01014                      uTotTdcHitToUseLeft, fVftxPar->uUseCoarseCorrectedTime )
01015                  - fVftxInputEvent->GetCalibratedTime(
01016                        fPar->uTdcPlasticLeftTime[uPlasticIndex],
01017                        fPar->uChannelPlasticLeftTime[uPlasticIndex],
01018                        iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01019                 )*fPar->dToTGainListLeft[uPlasticIndex]/1000.0
01020 //                - fPar->dOffsetListLeft[uPlasticIndex]/1000.0
01021                 - fPar->dTotOffsetListLeft[uPlasticIndex],
01022                (Double_t)uPlasticIndex );
01023 
01024             // Check last time info there => valid hit!
01025             if( kTRUE == fVftxInputEvent->IsHitThere(
01026                         fPar->uTdcPlasticRightTot[uPlasticIndex],
01027                         fPar->uChannelPlasticRightTot[uPlasticIndex],
01028                         iHitIndex ) )
01029             {
01030                // Check first if the main reference TDC is defined in VFTX parameter
01031                if( -1 < fVftxPar->iMainReferenceTdc )
01032                {
01033                   // Then Check if the 1st reference channel is defined for this TDC
01034                   if( -1 < fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc] )
01035                      if( kTRUE == fVftxInputEvent->IsHitThere(
01036                                     fVftxPar->iMainReferenceTdc,
01037                                     fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc],
01038                                     0 ) )
01039                      {
01040                         fReference1ProfilePlastics->Fill(
01041                         ( fVftxInputEvent->GetCalibratedTime(
01042                               fVftxPar->iMainReferenceTdc,
01043                               fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc],
01044                               0, fVftxPar->uUseCoarseCorrectedTime )
01045                             - ( fVftxInputEvent->GetCalibratedTime(
01046                                     fPar->uTdcPlasticRightTime[uPlasticIndex],
01047                                     fPar->uChannelPlasticRightTime[uPlasticIndex],
01048                                     iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01049                               + fVftxInputEvent->GetCalibratedTime(
01050                                     fPar->uTdcPlasticLeftTime[uPlasticIndex],
01051                                     fPar->uChannelPlasticLeftTime[uPlasticIndex],
01052                                     iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01053                               + fPar->dOffsetListRight[uPlasticIndex]
01054                               + fPar->dOffsetListLeft[uPlasticIndex]
01055                               - ( 1 == fVftxPar->uTdcOffsetEnable? dAutoOffsetSum : 0.0)
01056                             )/2.
01057                                 ),
01058                          (Double_t)uPlasticIndex );
01059                         // Then Check if the 2nd reference channel is defined for this TDC
01060                         if( -1 < fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc] )
01061                            if( kTRUE == fVftxInputEvent->IsHitThere(
01062                                           fVftxPar->iMainReferenceTdc,
01063                                           fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc],
01064                                           0 ) )
01065                            {
01066                               Double_t dRefMeanSum = 0.;
01067                               dRefMeanSum += fVftxInputEvent->GetCalibratedTime(
01068                                                 fVftxPar->iMainReferenceTdc,
01069                                                 fVftxPar->iVftxReference1Channel[fVftxPar->iMainReferenceTdc],
01070                                                 0, fVftxPar->uUseCoarseCorrectedTime );
01071                               dRefMeanSum += fVftxInputEvent->GetCalibratedTime(
01072                                                 fVftxPar->iMainReferenceTdc,
01073                                                 fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc],
01074                                                 0, fVftxPar->uUseCoarseCorrectedTime );
01075                               fMeanRefProfilePlastics->Fill(
01076                                ( ( dRefMeanSum/2. -
01077                                    ( fVftxInputEvent->GetCalibratedTime(
01078                                            fPar->uTdcPlasticRightTime[uPlasticIndex],
01079                                            fPar->uChannelPlasticRightTime[uPlasticIndex],
01080                                            iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01081                                      + fVftxInputEvent->GetCalibratedTime(
01082                                            fPar->uTdcPlasticLeftTime[uPlasticIndex],
01083                                            fPar->uChannelPlasticLeftTime[uPlasticIndex],
01084                                            iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01085                                      + fPar->dOffsetListRight[uPlasticIndex]
01086                                      + fPar->dOffsetListLeft[uPlasticIndex]
01087                                      - ( 1 == fVftxPar->uTdcOffsetEnable? dAutoOffsetSum : 0.0)
01088                                    )/2.
01089                                        ) ),
01090                                 (Double_t)uPlasticIndex );
01091                            }
01092                      }
01093                   // Then Check if the 2nd reference channel is defined for this TDC
01094                   if( -1 < fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc] )
01095                      if( kTRUE == fVftxInputEvent->IsHitThere(
01096                                     fVftxPar->iMainReferenceTdc,
01097                                     fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc],
01098                                     0 ) )
01099                      {
01100                         fReference2ProfilePlastics->Fill(
01101                          ( fVftxInputEvent->GetCalibratedTime(
01102                                fVftxPar->iMainReferenceTdc,
01103                                fVftxPar->iVftxReference2Channel[fVftxPar->iMainReferenceTdc],
01104                                0, fVftxPar->uUseCoarseCorrectedTime )
01105                            - ( fVftxInputEvent->GetCalibratedTime(
01106                                      fPar->uTdcPlasticRightTime[uPlasticIndex],
01107                                      fPar->uChannelPlasticRightTime[uPlasticIndex],
01108                                      iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01109                                + fVftxInputEvent->GetCalibratedTime(
01110                                      fPar->uTdcPlasticLeftTime[uPlasticIndex],
01111                                      fPar->uChannelPlasticLeftTime[uPlasticIndex],
01112                                      iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01113                                + fPar->dOffsetListRight[uPlasticIndex]
01114                                + fPar->dOffsetListLeft[uPlasticIndex]
01115                                - ( 1 == fVftxPar->uTdcOffsetEnable? dAutoOffsetSum : 0.0)
01116                              )/2.
01117                           ),
01118                           (Double_t)uPlasticIndex );
01119                      }
01120                } // if( -1 < fVftxPar->iMainReferenceTdc )
01121 
01122                // Store time including offset!
01123                hitCurrent.dTimeLeft  = fVftxInputEvent->GetCalibratedTime(
01124                                              fPar->uTdcPlasticLeftTime[uPlasticIndex],
01125                                              fPar->uChannelPlasticLeftTime[uPlasticIndex],
01126                                              iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01127                                        + fPar->dOffsetListLeft[uPlasticIndex]
01128                                        - dVftxOtherOffset
01129                                        - CLOCK_TIME*( 1 == fVftxPar->uTdcOffsetEnable?
01130                                              fVftxPar->iAutomaticTdcOffset[ fPar->uTdcPlasticLeftTime[ uPlasticIndex] ] : 0.0);
01131 
01132                // Store tot including offset and gain!
01133                // Missing the offset for the coarse counter overflow!!!!
01134                hitCurrent.dTotLeft   =
01135                   ( fVftxInputEvent->GetCalibratedTime(
01136                         fPar->uTdcPlasticLeftTot[uPlasticIndex],
01137                         fPar->uChannelPlasticLeftTot[uPlasticIndex],
01138                         uTotTdcHitToUseLeft, fVftxPar->uUseCoarseCorrectedTime )
01139                    - fVftxInputEvent->GetCalibratedTime(
01140                          fPar->uTdcPlasticLeftTime[uPlasticIndex],
01141                          fPar->uChannelPlasticLeftTime[uPlasticIndex],
01142                          iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01143                   )*fPar->dToTGainListLeft[uPlasticIndex]
01144 //                  - fPar->dOffsetListLeft[uPlasticIndex]
01145                   - fPar->dTotOffsetListLeft[uPlasticIndex]*1000.0;
01146                //Store time including offset!
01147                hitCurrent.dTimeRight = fVftxInputEvent->GetCalibratedTime(
01148                                              fPar->uTdcPlasticRightTime[uPlasticIndex],
01149                                              fPar->uChannelPlasticRightTime[uPlasticIndex],
01150                                              iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01151                                        + fPar->dOffsetListRight[uPlasticIndex]
01152                                        - dVftxOtherOffset
01153                                        - CLOCK_TIME*( 1 == fVftxPar->uTdcOffsetEnable?
01154                                              fVftxPar->iAutomaticTdcOffset[ fPar->uTdcPlasticRightTime[uPlasticIndex] ] : 0.0);
01155                // Store tot including offset and gain!
01156                // Missing the offset for the coarse counter overflow!!!!
01157                hitCurrent.dTotRight  =
01158                   ( fVftxInputEvent->GetCalibratedTime(
01159                          fPar->uTdcPlasticRightTot[uPlasticIndex],
01160                          fPar->uChannelPlasticRightTot[uPlasticIndex],
01161                          uTotTdcHitToUseRight, fVftxPar->uUseCoarseCorrectedTime )
01162                    - fVftxInputEvent->GetCalibratedTime(
01163                          fPar->uTdcPlasticRightTime[uPlasticIndex],
01164                          fPar->uChannelPlasticRightTime[uPlasticIndex],
01165                          iHitIndex, fVftxPar->uUseCoarseCorrectedTime )
01166                   )*fPar->dToTGainListRight[uPlasticIndex]
01167 //                  - fPar->dOffsetListRight[uPlasticIndex]
01168                   - fPar->dTotOffsetListRight[uPlasticIndex]*1000.0;
01169 
01170                ((fOutputEvent->fEvents[0]).fHits[uPlasticIndex]).push_back( hitCurrent );
01171 
01172                hitCurrent.Clear();
01173             }// if valid trailing edge also on right side => valid hit!
01174          } // if valid trailing edge on left side
01175       } // if( valid on both ends!)
01176 
01177       // Too close Hit merging attempt
01178       // For now go on only if both were merged
01179       // If only one => stop building hits for this event
01180       if( iHitIndex + 1 == (Int_t)uTotTdcHitToUseLeft )
01181       {
01182          if( iHitIndex + 1 == (Int_t)uTotTdcHitToUseRight )
01183             iHitIndex = uTotTdcHitToUseLeft;
01184             else break;
01185       }
01186       else if( iHitIndex + 1 == (Int_t)uTotTdcHitToUseRight )
01187          break;
01188    } // for TdcDataIndex in both end
01189 
01190    return;
01191 }
01192 /**************************************************************************************************/
01193 void TPlasticsProc::ProcessGet4v10PlasticSingle( UInt_t uEventIndex , UInt_t uPlasticIndex, Double_t dGet4v10OtherOffset )
01194 {
01195    // Test if both time are there
01196    for( UInt_t uHitIndex = 0;
01197          uHitIndex < ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01198                          fGet4Boards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).
01199                            fHits[ fPar->uChannelPlasticLeftTime[uPlasticIndex] ].size() );
01200          uHitIndex++)
01201    {
01202       // Left tot always there as GET4 unpack define hit with time+tot!
01203       fTotLeftPlastics->Fill(
01204           ( ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01205                    fGet4Boards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).
01206                      fHits[ fPar->uChannelPlasticLeftTime[uPlasticIndex] ] )[uHitIndex].GetTot()
01207           )*fPar->dToTGainListLeft[uPlasticIndex]/1000.0
01208           - fPar->dTotOffsetListLeft[uPlasticIndex],
01209          (Double_t)uPlasticIndex );
01210 
01211       // Store time including offset!
01212       hitCurrent.dTimeLeft  = ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01213                                  fGet4Boards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).
01214                                    fHits[ fPar->uChannelPlasticLeftTime[uPlasticIndex] ] )[uHitIndex].GetHitTime()
01215                               + fPar->dOffsetListLeft[uPlasticIndex]
01216                               - dGet4v10OtherOffset;
01217 
01218       // Store tot including offset and gain!
01219       // Missing the offset for the coarse counter overflow!!!!
01220       hitCurrent.dTotLeft   =
01221          ( ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01222                   fGet4Boards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).
01223                     fHits[ fPar->uChannelPlasticLeftTime[uPlasticIndex] ] )[uHitIndex].GetTot()
01224          )*fPar->dToTGainListLeft[uPlasticIndex]
01225          - fPar->dTotOffsetListLeft[uPlasticIndex]*1000.0;
01226 
01227       ((fOutputEvent->fEvents[uEventIndex]).fHits[uPlasticIndex]).push_back( hitCurrent );
01228 
01229       hitCurrent.Clear();
01230    } // for TdcDataIndex in both end
01231 
01232    return;
01233 }
01234 void TPlasticsProc::ProcessGet4v10PlasticDouble( UInt_t uEventIndex , UInt_t uPlasticIndex, Double_t dGet4v10OtherOffset )
01235 {
01236    // Test if both time are there
01237    for( UInt_t uHitIndex = 0;
01238          uHitIndex < ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01239                          fGet4Boards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).
01240                            fHits[ fPar->uChannelPlasticLeftTime[uPlasticIndex] ].size() ) &&
01241          uHitIndex < ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01242                          fGet4Boards[ fPar->uTdcPlasticRightTime[uPlasticIndex] ] ).
01243                            fHits[ fPar->uChannelPlasticRightTime[uPlasticIndex] ].size() );
01244          uHitIndex++)
01245    {
01246       fBeamProfilePlasticsTime->Fill(
01247                             ( ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01248                                     fGet4Boards[ fPar->uTdcPlasticRightTime[uPlasticIndex] ] ).
01249                                       fHits[ fPar->uChannelPlasticRightTime[uPlasticIndex] ] )[uHitIndex].GetHitTime()
01250                              -( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01251                                     fGet4Boards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).
01252                                       fHits[ fPar->uChannelPlasticLeftTime[uPlasticIndex] ] )[uHitIndex].GetHitTime()
01253                              + fPar->dOffsetListRight[uPlasticIndex]
01254                              - fPar->dOffsetListLeft[uPlasticIndex]
01255                                   ),
01256             ((Double_t)uPlasticIndex) +0.00001 );
01257 
01258       // Right Tot always there as GET4 unpack define hit with time+tot!
01259       fTotRightPlastics->Fill(
01260           ( ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01261                    fGet4Boards[ fPar->uTdcPlasticRightTime[uPlasticIndex] ] ).
01262                      fHits[ fPar->uChannelPlasticRightTime[uPlasticIndex] ] )[uHitIndex].GetTot()
01263           )*fPar->dToTGainListRight[uPlasticIndex]/1000.0
01264           - fPar->dTotOffsetListRight[uPlasticIndex],
01265          (Double_t)uPlasticIndex );
01266       // Left tot always there as GET4 unpack define hit with time+tot!
01267       fTotLeftPlastics->Fill(
01268           ( ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01269                    fGet4Boards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).
01270                      fHits[ fPar->uChannelPlasticLeftTime[uPlasticIndex] ] )[uHitIndex].GetTot()
01271           )*fPar->dToTGainListLeft[uPlasticIndex]/1000.0
01272           - fPar->dTotOffsetListLeft[uPlasticIndex],
01273          (Double_t)uPlasticIndex );
01274 
01275       // Store time including offset!
01276       hitCurrent.dTimeLeft  = ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01277                                  fGet4Boards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).
01278                                    fHits[ fPar->uChannelPlasticLeftTime[uPlasticIndex] ] )[uHitIndex].GetHitTime()
01279                               + fPar->dOffsetListLeft[uPlasticIndex]
01280                               - dGet4v10OtherOffset;
01281 
01282       // Store tot including offset and gain!
01283       // Missing the offset for the coarse counter overflow!!!!
01284       hitCurrent.dTotLeft   =
01285          ( ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01286                   fGet4Boards[ fPar->uTdcPlasticLeftTime[uPlasticIndex] ] ).
01287                     fHits[ fPar->uChannelPlasticLeftTime[uPlasticIndex] ] )[uHitIndex].GetTot()
01288          )*fPar->dToTGainListLeft[uPlasticIndex]
01289          - fPar->dTotOffsetListLeft[uPlasticIndex]*1000.0;
01290       //Store time including offset!
01291       hitCurrent.dTimeRight = ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01292                                  fGet4Boards[ fPar->uTdcPlasticRightTime[uPlasticIndex] ] ).
01293                                    fHits[ fPar->uChannelPlasticRightTime[uPlasticIndex] ] )[uHitIndex].GetHitTime()
01294                               + fPar->dOffsetListRight[uPlasticIndex]
01295                               - dGet4v10OtherOffset;
01296       // Store tot including offset and gain!
01297       // Missing the offset for the coarse counter overflow!!!!
01298       hitCurrent.dTotRight  =
01299          ( ( ( (fGet4v1InputEvent->fEvents[uEventIndex]).
01300                   fGet4Boards[ fPar->uTdcPlasticRightTime[uPlasticIndex] ] ).
01301                     fHits[ fPar->uChannelPlasticRightTime[uPlasticIndex] ] )[uHitIndex].GetTot()
01302          )*fPar->dToTGainListRight[uPlasticIndex]
01303          - fPar->dTotOffsetListRight[uPlasticIndex]*1000.0;
01304 
01305       ((fOutputEvent->fEvents[uEventIndex]).fHits[uPlasticIndex]).push_back( hitCurrent );
01306 
01307       hitCurrent.Clear();
01308    } // for TdcDataIndex in both end
01309 
01310    return;
01311 }
01312 /**************************************************************************************************/

Generated on Tue Dec 10 2013 04:52:20 for ROCsoft by  doxygen 1.7.1