SpillTimeFinder Class Reference

#include <SpillTimeFinder.h>

Inheritance diagram for SpillTimeFinder:
CfgPromptConfigurable

List of all members.

Classes

struct  Cleaner

Public Member Functions

Bool_t DataIsAvailable (const VldContext &context)
const SpillTimeNDGetRecentSpill (const VldContext &context)
const SpillTimeNDGetNextSpill (const VldContext &context)
const SpillTimeNDGetNearestSpill (const VldContext &context)
VldTimeStamp GetTimeOfRecentSpill (const VldContext &cx)
VldTimeStamp GetTimeOfNextSpill (const VldContext &cx)
VldTimeStamp GetTimeOfNearestSpill (const VldContext &cx)
Double_t GetTimeToRecentSpill (const VldContext &cx)
Double_t GetTimeToNextSpill (const VldContext &cx)
Double_t GetTimeToNearestSpill (const VldContext &cx)
VldTimeStamp GetLastQueryRangeBeg ()
VldTimeStamp GetLastQueryRangeEnd ()
 SpillTimeFinder ()
virtual ~SpillTimeFinder ()

Static Public Member Functions

static SpillTimeFinderInstance ()
static void SwitchOffsetsOff ()
static void SwitchOffsetsOn ()

Private Member Functions

virtual void ConfigModified ()
void FindClosestEntries (const VldContext &inContext, const SpillTimeND *&outPrev, const SpillTimeND *&outNext)
void FindBestRows (DbiResultPtr< SpillTimeND > &table, const VldTimeStamp &time, const SpillTimeND *&outPrev, const SpillTimeND *&outNext)
Double_t GetOffset (const VldContext &context, Int_t task)
Double_t GetOffsetNDNuToFDNu (const VldContext &context)
Double_t GetOffsetKickerToNDNu (const VldContext &context)
Double_t GetOffsetSgateToNDNu (const VldContext &context)
Double_t GetTimeDifference (const VldTimeStamp &t1, const VldTimeStamp &v2)
void ResetCache ()

Private Attributes

DbiResultPtr< SpillTimeNDfCurrentTable
DbiResultPtr< SpillTimeND > * fPrevTable
DbiResultPtr< SpillTimeND > * fNextTable
DbiResultPtr
< SpillTimeCalibration
fCalibrationNDNuToFDNu
DbiResultPtr
< SpillTimeCalibration
fCalibrationKickerToNDNu
DbiResultPtr
< SpillTimeCalibration
fCalibrationSgateToNDNu
SpillTimeND fLastQueryPrev
SpillTimeND fLastQueryNext
VldTimeStamp fLastQueryRangeBeg
VldTimeStamp fLastQueryRangeEnd
Int_t fTask

Static Private Attributes

static Bool_t fSwitchOffsetsOff = false
static SpillTimeFinderfgInstance = 0

Friends

struct Cleaner

Detailed Description

Definition at line 10 of file SpillTimeFinder.h.


Constructor & Destructor Documentation

SpillTimeFinder::SpillTimeFinder (  ) 

Definition at line 37 of file SpillTimeFinder.cxx.

References CfgPromptConfigurable::InitializeConfig(), and Registry::Set().

Referenced by Instance().

00037                                  :
00038   fPrevTable(0),
00039   fNextTable(0),
00040   fLastQueryPrev(kSpillTime_Future),
00041   fLastQueryNext(kSpillTime_VeryOld)
00042 {
00043   Registry r;
00044   r.Set("Task",3);
00045   InitializeConfig(r);
00046 }

SpillTimeFinder::~SpillTimeFinder (  )  [virtual]

Definition at line 48 of file SpillTimeFinder.cxx.

References fNextTable, and fPrevTable.

00049 {
00050   if(fPrevTable) delete fPrevTable;
00051   if(fNextTable) delete fNextTable;
00052 }


Member Function Documentation

void SpillTimeFinder::ConfigModified ( void   )  [private, virtual]

Implements CfgPromptConfigurable.

Definition at line 60 of file SpillTimeFinder.cxx.

References fTask, Registry::Get(), CfgPromptConfigurable::GetConfig(), Msg::kWarning, MSG, and ResetCache().

00061 {
00062   bool ok = true;
00063   ok = ok && GetConfig().Get("Task",fTask);
00064   if (!ok) {
00065     MSG("SpillTimeFinder",Msg::kWarning)
00066       << "Problem configuring SpillTimeFinder." << endl;  
00067   }
00068   ResetCache();
00069 }

Bool_t SpillTimeFinder::DataIsAvailable ( const VldContext context  ) 

Return true if there is a database table that covers this period of time. Return false if this period is not covered.

Definition at line 71 of file SpillTimeFinder.cxx.

References VldTimeStamp::Add(), fCurrentTable, fTask, GetOffset(), VldContext::GetSimFlag(), VldContext::GetTimeStamp(), DbiResultPtr< T >::GetValidityRec(), DbiValidityRec::IsGap(), Detector::kNear, and DbiResultPtr< T >::NewQuery().

Referenced by NtpMaker::FillSpillInfo(), and DataQualityInterface::ProcessBeamStatus().

00072 {
00077 
00078   // This is the relevant time at the kicker fire $74 MIBS
00079   VldTimeStamp kickerTime = context.GetTimeStamp();
00080   Double_t offset = GetOffset(context,fTask);
00081   kickerTime.Add(-offset); // Move requested time back to get kicker time.
00082 
00083   VldContext kickerContext(Detector::kNear,
00084                            context.GetSimFlag(),
00085                            kickerTime);
00086 
00087   // Query on the current context.
00088   fCurrentTable.NewQuery(kickerContext,fTask);
00089   
00090   // <Refresh cache here if any.>
00091 
00092 
00093   const DbiValidityRec* myVrec = fCurrentTable.GetValidityRec();
00094   if ( myVrec->IsGap () ) return false; // No data available
00095   
00096   return true;  // Data available.
00097 }

void SpillTimeFinder::FindBestRows ( DbiResultPtr< SpillTimeND > &  table,
const VldTimeStamp time,
const SpillTimeND *&  outPrev,
const SpillTimeND *&  outNext 
) [private]

Find the best match within the confines of the current table. Return the two entries corresponding to time>=outPrev, time<outNext or NULL if the entry doesn't exist in the table.

Definition at line 621 of file SpillTimeFinder.cxx.

References DbiResultPtr< T >::GetNumRows(), DbiResultPtr< T >::GetRow(), GetTimeDifference(), and SpillTimeND::GetTimeStamp().

Referenced by FindClosestEntries().

00625 {
00631   
00632   //
00633   outPrev = outNext = 0;
00634 
00635   UInt_t low = 0;
00636   UInt_t high= table.GetNumRows();
00637 
00638   if(high <= 0) return; // No rows available. (Query must have failed.)
00639 
00640   // Can't depend on table being sorted, so I'll just brute-force the mother.
00641 
00642   double dt;
00643   double minneg=1e99;
00644   double minpos=1e99;
00645 
00646 
00647   for(UInt_t i=low;i<high;i++) {
00648     const SpillTimeND* row = table.GetRow(i);
00649     dt = GetTimeDifference(row->GetTimeStamp(),time);
00650     if(dt <=0) {
00651       if( (-dt)<minneg ) {
00652         minneg = -dt;
00653         outPrev = row;
00654       }
00655     } else {
00656       if( (dt)<minpos ) {
00657         minpos = -dt;
00658         outNext = row;
00659       }
00660     }
00661   }
00662 
00663   // Throw out any solutions that are just plain stupid.
00664   const VldTimeStamp tooEarly(1975,0,0,0,0,0,0);
00665 
00666   if(outPrev) 
00667         if(outPrev->GetTimeStamp() < tooEarly) outPrev = 0;
00668   if(outNext)
00669         if(outNext->GetTimeStamp() < tooEarly) outNext = 0;
00670 }

void SpillTimeFinder::FindClosestEntries ( const VldContext inContext,
const SpillTimeND *&  outPrev,
const SpillTimeND *&  outNext 
) [private]

Search the available databases for spill times close to the given context. Return spill records for time <= context, and > context.

Definition at line 101 of file SpillTimeFinder.cxx.

References VldTimeStamp::Add(), VldContext::AsString(), VldTimeStamp::AsString(), fCurrentTable, FindBestRows(), fLastQueryNext, fLastQueryPrev, fLastQueryRangeBeg, fLastQueryRangeEnd, fNextTable, Form(), fPrevTable, fTask, DbiResultPtr< T >::GetNumRows(), GetOffset(), VldTimeStamp::GetSec(), VldContext::GetSimFlag(), VldRange::GetTimeEnd(), SpillTimeND::GetTimeStamp(), VldContext::GetTimeStamp(), VldRange::GetTimeStart(), DbiResultPtr< T >::GetValidityRec(), DbiValidityRec::GetVldRange(), Dbi::kAnyTask, Msg::kDebug, Msg::kError, SimFlag::kMC, Detector::kNear, kSpillTime_VeryOld, Msg::kVerbose, MSG, n, and DbiResultPtr< T >::NewQuery().

Referenced by GetNearestSpill(), GetNextSpill(), and GetRecentSpill().

00104 {
00105   
00110 
00111   if(context.GetSimFlag()==SimFlag::kMC) {
00112     MSG("SpillTimeFinder",Msg::kError)
00113       << "MC context in SpillTimeFinder request. This will all end in tears, I know it..." << endl;
00114   }
00115 
00116   outPrev = outNext = 0;
00117 
00118   // This is the relevant time at the kicker fire $74 MIBS
00119   VldTimeStamp kickerTime = context.GetTimeStamp();
00120   Double_t offset = GetOffset(context,fTask);
00121   kickerTime.Add(-offset); // Move requested time back to get kicker time.
00122 
00123   VldContext kickerContext(Detector::kNear,
00124                            context.GetSimFlag(),
00125                            kickerTime );
00126 
00127   MSG("SpillTimeFinder",Msg::kVerbose) 
00128     << "STF Query: " << context.AsString() << endl;
00129   MSG("SpillTimeFinder",Msg::kVerbose) 
00130     << "Kicker Context: " << kickerContext.AsString() << endl;
00131   
00132   
00133   // Optimization:
00134   // If this request is in (prev,next) of the last call,
00135   // the same prev,next results should still apply.
00136   if(kickerTime < fLastQueryNext.GetTimeStamp()) {
00137     if(kickerTime > fLastQueryPrev.GetTimeStamp()) {
00138       // Match! Just give the last results.
00139       MSG("SpillTimeFinder",Msg::kDebug) << "Using last call's cached result." << endl;
00140       outPrev = &fLastQueryPrev;
00141       outNext = &fLastQueryNext;
00142       MSG("SpillTimeFinder",Msg::kVerbose) << "Result: " << outPrev->GetTimeStamp().AsString() 
00143                                            << "  " << outNext->GetTimeStamp().AsString() << endl;
00144 
00145       return;
00146     }     
00147   }
00148 
00149 
00150   // Query on the current context.
00151   fCurrentTable.NewQuery(kickerContext,fTask);
00152 
00153   // <Refresh cache here if any.>
00154 
00155   MSG("SpillTimeFinder",Msg::kDebug) << "Requested context, got " << fCurrentTable.GetNumRows() << " rows" << endl;  
00156   // See what rows we can get.
00157   FindBestRows(fCurrentTable,kickerTime, outPrev, outNext);
00158 
00159   MSG("SpillTimeFinder",Msg::kVerbose) << "First query: Prev = "
00160                                        << ((outPrev==0)?"<null>":outPrev->GetTimeStamp().AsString())
00161                                        << "   " 
00162                                        << ((outNext==0)?"<null>":outNext->GetTimeStamp().AsString())
00163                                        << endl;
00164 
00165 
00166 
00167   // Are we done?
00168   if((outPrev)&&(outNext)) return;
00169 
00170 
00172   // Ok, we need to hunt for more info.
00173 
00174   // Default: if we didn't find a table at all, set the search param
00175   // as the requested param.
00176   VldTimeStamp curStartVldTime;
00177   VldTimeStamp curStopVldTime;
00178 
00179   // Get the time limits of the current query.
00180   const DbiValidityRec* curValidity = fCurrentTable.GetValidityRec();
00181   if(curValidity) {
00182     // Use the start and stop times here.
00183     curStartVldTime = curValidity->GetVldRange().GetTimeStart();
00184     curStopVldTime = curValidity->GetVldRange().GetTimeEnd();
00185   } else {
00186     curStartVldTime = kickerTime;
00187     curStopVldTime = kickerTime;
00188   }
00189 
00190   fLastQueryRangeBeg = curStartVldTime; 
00191   fLastQueryRangeBeg.Add(GetOffset(context,fTask));
00192   fLastQueryRangeEnd = curStopVldTime;  
00193   fLastQueryRangeEnd.Add(GetOffset(context,fTask));
00194 
00195 
00196   // Did we miss the early side?
00197   while(outPrev ==0) {
00198     if(curStartVldTime < VldTimeStamp(1975,0,0,0,0,0)){
00199       MSG("SpillTimeFinder",Msg::kDebug) << "Hit start of vld range. Hit starting bookend (<1975 AD)" << endl;
00200       outPrev = &kSpillTime_VeryOld;
00201       break;
00202     }
00203     
00204     MSG("SpillTimeFinder",Msg::kDebug) << "Hit start of vld range.. looking earlier..." << endl;
00205     // Look for an earlier table.
00206     const char* sqltxt = Form("(TIMEEND<='%s') and (TASK=%d) and (DETECTORMASK & %d) and (SIMMASK & %d) order by TIMEEND desc limit 1",
00207                               curStartVldTime.AsString("s"),
00208                               fTask,
00209                               kickerContext.GetDetector(),
00210                               kickerContext.GetSimFlag() );
00211     MSG("SpillTimeFinder",Msg::kDebug) << "Request: " << sqltxt << endl;
00212     
00213     DbiSqlContext myContext(sqltxt);
00214     if(fPrevTable) {
00215       fPrevTable->NewQuery(myContext,Dbi::kAnyTask);
00216     } else {
00217       fPrevTable = new DbiResultPtr<SpillTimeND>("SPILLTIMEND",myContext);
00218     }
00219     
00220     if(fPrevTable->GetNumRows()>0) {
00221       // Found a table with data! Look for the best row..
00222       MSG("SpillTimeFinder",Msg::kDebug)
00223         << "...earlier table gotten. Rows: "
00224         << fPrevTable->GetNumRows() << endl;
00225       const SpillTimeND* n;    
00226       FindBestRows(*fPrevTable,kickerTime,outPrev,n);
00227       if ( n > 0 ) {
00228         MSG("SpillTimeFinder",Msg::kError)
00229           << "Conflict! Found 'next' time in 'previous' table: "
00230           << n->GetTimeStamp().AsString() << endl;    
00231       }
00232       
00233     } else {
00234       VldTimeStamp nexttry = fPrevTable->GetValidityRec()->GetVldRange().GetTimeStart();
00235       if(nexttry >= curStartVldTime) nexttry = VldTimeStamp(curStartVldTime.GetSec()-1,0); // protect against infinite loop.
00236       curStartVldTime = nexttry;
00237     }
00238   }
00239  
00240 
00241   // Did we miss the late side?
00242   while(outNext ==0) {
00243     MSG("SpillTimeFinder",Msg::kDebug) << "Hit end of vld range.. looking later..." << endl;
00244     const char* sqltxt = Form("(TIMESTART>='%s') and (TASK=%d) and (DETECTORMASK & %d) and (SIMMASK & %d) order by TIMESTART asc limit 1",
00245                               curStopVldTime.AsString("s"),
00246                               fTask,
00247                               kickerContext.GetDetector(),
00248                               kickerContext.GetSimFlag()
00249                               );
00250     
00251     MSG("SpillTimeFinder",Msg::kDebug) << "Request: " << sqltxt << endl;
00252     
00253     DbiSqlContext myContext(sqltxt);
00254     if(fNextTable)
00255       fNextTable->NewQuery(myContext,Dbi::kAnyTask);
00256     else 
00257       fNextTable = new DbiResultPtr<SpillTimeND>("SPILLTIMEND",myContext);
00258     
00259     if(fNextTable->GetNumRows()>0) {
00260       // Got a table with data! Look for the best row:
00261       MSG("SpillTimeFinder",Msg::kDebug) << "...later table gotten. Rows: " << fNextTable->GetNumRows() << endl;
00262       const SpillTimeND* p;
00263       
00264       FindBestRows(*fNextTable,kickerTime,p,outNext);
00265       if ( p > 0 ) {
00266         MSG("SpillTimeFinder",Msg::kError)
00267           << "Conflict! Found 'previous' time in 'next' table: "
00268           << p->GetTimeStamp().AsString() << endl;
00269       }
00270       
00271     } else {
00272      VldTimeStamp nexttry = fNextTable->GetValidityRec()->GetVldRange().GetTimeStart();
00273      if(nexttry <= curStopVldTime) nexttry = VldTimeStamp(curStopVldTime.GetSec()+1,0); // protect against infinite loop.
00274      curStopVldTime = nexttry;
00275     }
00276   }
00277 
00278   // Sanity check.
00279   if ( outPrev->GetTimeStamp() > kickerTime ) {
00280     MSG("SpillTimeFinder",Msg::kError)
00281       << "Sanity check failed. Prev time is late!";
00282   }
00283   if ( outNext->GetTimeStamp() < kickerTime ) {
00284     MSG("SpillTimeFinder",Msg::kError)
00285       << "Sanity check failed. Next time is early!";
00286   }
00287 
00288   // Load the cache.
00289   fLastQueryPrev = *outPrev;
00290   fLastQueryNext = *outNext;
00291 
00292   MSG("SpillTimeFinder",Msg::kVerbose) << "Result: " << outPrev->GetTimeStamp().AsString() 
00293                                        << "  " << outNext->GetTimeStamp().AsString() << endl;
00294 
00295 }

VldTimeStamp SpillTimeFinder::GetLastQueryRangeBeg (  )  [inline]

Definition at line 36 of file SpillTimeFinder.h.

References fLastQueryRangeBeg.

00036 { return fLastQueryRangeBeg; };

VldTimeStamp SpillTimeFinder::GetLastQueryRangeEnd (  )  [inline]

Definition at line 37 of file SpillTimeFinder.h.

References fLastQueryRangeEnd.

Referenced by SpillTimeCreateKeyFile().

00037 { return fLastQueryRangeEnd; };

const SpillTimeND & SpillTimeFinder::GetNearestSpill ( const VldContext context  ) 

Definition at line 708 of file SpillTimeFinder.cxx.

References VldTimeStamp::Add(), FindClosestEntries(), fTask, GetOffset(), GetTimeDifference(), SpillTimeND::GetTimeStamp(), VldContext::GetTimeStamp(), and kSpillTime_VeryOld.

Referenced by MadDpAnalysis::CreatePAN(), MadTestAnalysis::CreatePAN(), MadPIDAnalysis::CreatePAN(), NtpMaker::FillSpillInfo(), GetTimeOfNearestSpill(), DataQualityInterface::ProcessBeamStatus(), and DataUtil::QueryBeamDB().

00709 {
00710   // return referenct to prev or next spill record, whichever is closer.
00711 
00712   const SpillTimeND* prev;
00713   const SpillTimeND* next;
00714   FindClosestEntries(context,prev,next);
00715 
00716   if((next==0)  && (prev==0)) return kSpillTime_VeryOld;
00717 
00718   if(!prev) return *next;
00719   if(!next) return *prev;
00720 
00721   double offset = GetOffset(context,fTask);
00722   VldTimeStamp prevts = prev->GetTimeStamp();
00723   prevts.Add(offset);
00724   VldTimeStamp nextts = next->GetTimeStamp();
00725   nextts.Add(offset);
00726 
00727   double dt1 = GetTimeDifference(prevts,context.GetTimeStamp());
00728   double dt2 = GetTimeDifference(nextts,context.GetTimeStamp());
00729 
00730   if(fabs(dt1)<=fabs(dt2)) {
00731     return *prev;
00732   }
00733   // else
00734   return *next; 
00735 }

const SpillTimeND & SpillTimeFinder::GetNextSpill ( const VldContext context  ) 

Definition at line 697 of file SpillTimeFinder.cxx.

References FindClosestEntries(), and kSpillTime_Future.

Referenced by GetTimeOfNextSpill().

00698 {
00699   // Return reference to next spill record.
00700 
00701   const SpillTimeND* prev;
00702   const SpillTimeND* next;
00703   FindClosestEntries(context,prev,next);
00704   if(next) return *next;
00705   else     return kSpillTime_Future;
00706 }

Double_t SpillTimeFinder::GetOffset ( const VldContext context,
Int_t  task 
) [private]

Compute the offset between the time in the database and the time the neutrinos appear to arrive at the ND or FD (depending on the context) For a complete description of all of this, see NuMI-NOTE-ANA-1090

Definition at line 299 of file SpillTimeFinder.cxx.

References fSwitchOffsetsOff, VldContext::GetDetector(), GetOffsetKickerToNDNu(), GetOffsetNDNuToFDNu(), GetOffsetSgateToNDNu(), Detector::kFar, Msg::kInfo, Detector::kNear, SpillTimeND::kTask_Vtm, and MAXMSG.

Referenced by DataIsAvailable(), FindClosestEntries(), GetNearestSpill(), GetTimeOfNearestSpill(), GetTimeOfNextSpill(), and GetTimeOfRecentSpill().

00300 {
00307 
00308   Double_t offset = 0;
00309 
00310   if (fSwitchOffsetsOff){
00311     if (Detector::kFar == context.GetDetector()){
00312       offset = 0.00244935653818216;
00313       MAXMSG("SpillTimeFinder",Msg::kInfo,20)
00314         << "SpillTime offsets turned off (FD)" << endl;
00315     }
00316     else{
00317       offset = 0.0;
00318       MAXMSG("SpillTimeFinder",Msg::kInfo,20)
00319         << "SpillTime offsets turned off (ND)" << endl;
00320     }
00321     return offset;
00322   }
00323   
00324   switch(context.GetDetector()) {
00325   case Detector::kFar:
00326     offset += GetOffsetNDNuToFDNu(context); //fallthrough
00327 
00328   case Detector::kNear:
00329     if(task == SpillTimeND::kTask_Vtm) {
00330       offset += GetOffsetSgateToNDNu(context);
00331     } else {
00332       offset += GetOffsetKickerToNDNu(context); // Get offset from kicker to ND
00333     }
00334     break;
00335   default:
00336     break;    
00337   }
00338   return offset;
00339 }

Double_t SpillTimeFinder::GetOffsetKickerToNDNu ( const VldContext context  )  [private]

Return time in seconds between the time the time the kicker gets MIBS $74, and the time that neutrinos hit the ND. which is what is recorded in the db

Constants are correct for late May, 2005.

Definition at line 513 of file SpillTimeFinder.cxx.

References fCalibrationKickerToNDNu, Form(), SpillTimeCalibration::GetDescription(), SpillTimeCalibration::GetError(), DbiResultPtr< T >::GetNumRows(), SpillTimeCalibration::GetOffset(), DbiResultPtr< T >::GetRow(), Msg::kInfo, Munits::microsecond, MSG, and DbiResultPtr< T >::NewQuery().

Referenced by GetOffset().

00514 {
00519 
00520   double offset = 0;
00521   fCalibrationKickerToNDNu.NewQuery(cx,1);
00522   if(fCalibrationKickerToNDNu.GetNumRows()>0) {
00523     static int sCall = 0;
00524     if ( sCall == 0 ) {
00525       MSG("SpillTimeFinder",Msg::kInfo)
00526         << "Initial Spill Timing Calibration constants (KickerToNDNu):" << endl;
00527     }
00528     for (unsigned int i = 0; i<fCalibrationKickerToNDNu.GetNumRows(); i++) {
00529       if ( sCall == 0 ) {
00530         MSG("SpillTimeFinder",Msg::kInfo) 
00531           << Form("%10.5e",fCalibrationKickerToNDNu.GetRow(i)->GetOffset())
00532           << " +- "
00533           << fCalibrationKickerToNDNu.GetRow(i)->GetError() 
00534           << "\t"
00535           << fCalibrationKickerToNDNu.GetRow(i)->GetDescription() << endl;
00536       }
00537       offset += fCalibrationKickerToNDNu.GetRow(i)->GetOffset();
00538     }
00539     sCall++;
00540     return offset;
00541   }
00542 
00543 
00545   const double kOffset = 
00546     // The following two values are taken directly from the 
00547     // TimeMaster control settings: SgateDelayRevs and SgateDelayFine:
00548     +19*11.064*Munits::microsecond // $74 MIBS until actual kicker fire, in MI turns
00549     +6.2*Munits::microsecond       // Fine tuning of SGATE window to $74
00550 
00551     // This tuning gets the first neutrinos at 1.5us into the SGATE
00552     // so add this on.  Error is +-0.1 us
00553     //+1.5*Munits::microsecond
00554     // Fixed 7/6/2006: retune based upon TOF study; move by -2.4 us
00555     -0.9*Munits::microsecond
00556 
00557     // This delay is not precise, because it does not truely wait an integer
00558     // number of turns, but rather waits for the NEXT n turns. Therefore, this 
00559     // delay is actually shorter. Tom Fitzpatrick measured this as being ~1us. (email 23/05/05)
00560     // Error not known. Estimate worst-case as 100%, or 1.0us
00561     - 1.0*Munits::microsecond
00562 
00563     // Now we have time between $74 and SGATE at the MCC. 
00564     // But, it takes finite time to propagate SGATE to the front end.
00565     // Tom Fitzpatrick measured this too (email 23/05/05) as being 450ns
00566     // Error is negligable compared to other offsets
00567     - 0.450*Munits::microsecond
00568 
00569   ;
00570 
00571   return kOffset;
00572 }

Double_t SpillTimeFinder::GetOffsetNDNuToFDNu ( const VldContext context  )  [private]

Return time in seconds between the time the spill hits the FD and the time it hits the ND (FD-ND)

The offsets in this code are best documented in the Time Of Flight

------------------------------------- FD latencies: --------------------------- /////////////////////// [7] //////////////////////////

Definition at line 341 of file SpillTimeFinder.cxx.

References fCalibrationNDNuToFDNu, Form(), SpillTimeCalibration::GetDescription(), SpillTimeCalibration::GetError(), DbiResultPtr< T >::GetNumRows(), SpillTimeCalibration::GetOffset(), DbiResultPtr< T >::GetRow(), VldContext::GetTimeStamp(), Msg::kInfo, Msg::kWarning, MAXMSG, Munits::microsecond, Munits::millisecond, MSG, Munits::nanosecond, and DbiResultPtr< T >::NewQuery().

Referenced by GetOffset().

00342 {
00348   // paper, DocDB 1870.
00349   
00350   // A note on signs:
00351   // A positive offset indicates the event will appear to be later at the FD
00352   // A negative offset indicates the event will appear to be earlier
00353   // Things that make the ND clock run slow will give a positive offset
00354   // Things that make the FD clock run slow will give a negative offset
00355 
00356   // New version - May 2008. Get offset calibration from DB, if available
00357   
00358   double offset = 0;
00359   fCalibrationNDNuToFDNu.NewQuery(cx,0);
00360   if(fCalibrationNDNuToFDNu.GetNumRows()>0) {
00361     static int sCall = 0;
00362     if ( sCall == 0 ) {
00363       MSG("SpillTimeFinder",Msg::kInfo)
00364         << "Initial Spill Timing Calibration constants (NDNuToFDNu):" << endl;
00365     }
00366     for (unsigned int i = 0; i<fCalibrationNDNuToFDNu.GetNumRows(); i++) {
00367       if ( sCall == 0 ) {
00368         MSG("SpillTimeFinder",Msg::kInfo) 
00369           << Form("%10.5e",fCalibrationNDNuToFDNu.GetRow(i)->GetOffset())
00370           << " +- "
00371           << fCalibrationNDNuToFDNu.GetRow(i)->GetError() 
00372           << "\t"
00373           << fCalibrationNDNuToFDNu.GetRow(i)->GetDescription() << endl;
00374       }
00375       offset += fCalibrationNDNuToFDNu.GetRow(i)->GetOffset();
00376     }
00377     sCall++;
00378     return offset;
00379   }
00380   
00381   MAXMSG("SpillTimeFinder",Msg::kWarning,1) << "Could not find SPILLTIMECALIBRATION table." << endl;
00382   MAXMSG("SpillTimeFinder",Msg::kWarning,1) << "Resorting to hard-coded calibration." << endl;
00383 
00384   const VldTimeStamp ndSwitch(2006,7,27,14,40,0,0); // Time Alfons changed ND cable.
00385   
00387   // ND-FD time of flight, 734,298.62 meters (Wes Smart), at c
00388   // Error is 0.7m, or 2.3ns.  
00389   offset += 2.44935653818215832*Munits::millisecond;
00390 
00391   //--------------------------------------
00392   // ND latencies
00393 
00395   // ND antenna:
00396   // OLD: total guess is 300ft, with same prop. speed as FD antenna.
00397   // This delays counter reset at ND
00398   // Error is +-60%, or 0.3us. (It's unlikely to be longer than 500ft, or shorter than 200)
00399   //+0.5*Munits::microsecond 
00400   // NEW April 5,2006. measurement by Chuck Andrews, reported by Andrew Rader
00401   // Fibre is in fact 765ft +/- 1% between the upstairs tranciever and U/G DAQ rack.
00402   // Andy gives the refrativce index as 1.4677 at 1300nm. 
00403   // Error is is 1%, or 12ns.
00404   // offset += 1.1415*Munits::microsecond;
00405   //
00406   // NEW Dec 6 2006:
00407   // Measured again with the BOLT. Got a new length
00408   // of 227.5m. Multiply by 1.470/3e8 s/m (the BOLT conversion constant).
00409   // So offset =  1.1155us +- 28ns
00410   // Take difference between measurements as the error: 
00411   //   1.1415-1.1155 = 0.027us
00412   offset += 1.1147*Munits::microsecond;
00413   
00415   // ND cable run from DAQ racks to GPS reciever:
00416   if(cx.GetTimeStamp() < ndSwitch) {
00417     // Cable run from ND tranciever to reciever (DAQ rack to timing rack)
00418     // Cat J sez: 5 segments of 32ns RG59 unioned together. +/-10% or 16 ns.
00419     offset += 0.160*Munits::microsecond;
00420   } else {
00421     // Alfons: we actually replaced 2m optical patch cable in the beams division crate
00422     // with a 50m optical patch cable, which leads now to the timing rack.
00423     // The chain of N x 10' RG58 is still in the cable tray, but not used.
00424     // Instead, the optical->elect converter is now in the timing rack and
00425     // connected with a (2+/-0.5) m RG58 to the GPS receiver    
00426     offset += 6.7*Munits::nanosecond // 2m electrical cable
00427       +0.2446*Munits::microsecond; // 50m optical cable
00428     // Another ~5ns of error here.
00429   }
00430   
00432   // ND timing system latency, MCC to MTM
00433   // But, it takes finite time to propagate CNTRST to the front end.
00434   // Tom Fitzpatrick measured this too (email 23/05/05) as being 450ns
00435   // Take reasonable error as +-50ns
00436   offset += 0.450*Munits::microsecond;
00437 
00439   // ND timing system latency, MINDER
00440   // The timing signals from the PMT to the timestamper is
00441   // finely timed in, so as to be less than 1 bucket. --Dave Reyna
00442   // Take offset as 10+-10ns
00443   offset += 10*Munits::nanosecond;
00444 
00446   // CandStrip reconstruction bias.
00447   // The ND measures the mean time of each strip, not the time of the first digit
00448   // The time of the first digit is a better estimate (compared to the FD), so correct for this
00449   // The bias varies as a function of pulse height, but mean bias is 7ns, +- 9ns rms
00450   offset += -7.0*Munits::nanosecond;
00451     
00454   // FD antenna: re-measured by Dave Saranen. Fibre length ~3415ft Email 23/05/05
00455   //
00456   // This delays counter reset at FD.
00457   // Error is on precision of the ping measurment: +-0.1 us
00458   // offset += -5.1*Munits::microsecond;
00459   
00460   // Much better measurement by Alec Habig, using our custom-bought BOLT unit.
00461   // Email Nov 8, 2006
00462   // Bolt gave a measurement of 1049m for the antenna run, and uses a group velocity of 
00463   // n_eff=1.470 and c=3e8
00464   // Disagreement between Orlando (FNAL) and BOLT measurement in for the
00465   // multimode fibre between VAX room and rescue room is 4m.  Error on BOLT is claimed to
00466   // be 2.5m.  Error on Orlando's measurement was 6m. 
00467   // Use differences in multimode fibre to see the error.
00468   offset += -5140*Munits::nanosecond;
00469   
00471   // FD timing system latency. Measured by using the RawVaTimingMonitorBlock
00472   // which gives a stable value of 639980381 TDC counts for the PPS fiducial tick.
00473   // Precision is +-1.5ns, but varies from run to run within ~100 ns.
00474   // Negative: the FD clock is even slower because of this latency.
00475   offset += -30.6547*Munits::microsecond;
00476   
00478   // New, Sept 2006:FD timing system fiducial cable length.
00479   // Makes the above measurement just a wee 
00480   // bit wrong. +-6ns (the rise time of the test pulse).
00481   // Negative.. it means the fiducial round-trip was even longer than we thought.
00482   offset += -153.0 * Munits::nanosecond; 
00483 
00485   // New Dec 2006: Average timing constants at the FD are not exactly
00486   // zero, but have an average offset of -27ns. This makes the events appear to be too early.
00487   // Error is half the offset, 14ns.
00488   offset += -27 * Munits::nanosecond;
00489   
00490   
00492   // Differences between ND and FD transceivers and recivers.
00493   // Difference = 0 +/- 20ns.
00494 
00495 
00497   // Detector effects. Total strip half-lenth plus WLS pigtail plus clear fibre readout
00498   // is 5m@ND and 8m@FD. (Very roughly; depends on the actual centroid of
00499   // accepted events vertices.) Thus, the events will appear at the FD about
00500   // 16ns later than events in the ND, due to the optical transit time.
00501   // Error on this is about 1m, or 5ns.
00502   offset += 16.0*Munits::nanosecond;
00503 
00504 
00506   // PMT transit time difference.
00507   // Both M16 and M64 have a 10.9ns transit time. Error is 1ns.
00508   
00509   
00510   return offset;
00511 }

Double_t SpillTimeFinder::GetOffsetSgateToNDNu ( const VldContext context  )  [private]

Get the time difference between when SGATE hits the front end and when the neutrinos hit the ND.

Definition at line 574 of file SpillTimeFinder.cxx.

References fCalibrationSgateToNDNu, Form(), SpillTimeCalibration::GetDescription(), SpillTimeCalibration::GetError(), DbiResultPtr< T >::GetNumRows(), SpillTimeCalibration::GetOffset(), DbiResultPtr< T >::GetRow(), VldTimeStamp::GetSec(), VldContext::GetTimeStamp(), Msg::kInfo, Munits::microsecond, MSG, and DbiResultPtr< T >::NewQuery().

Referenced by GetOffset().

00575 {
00580 
00581   double offset = 0;
00582   fCalibrationSgateToNDNu.NewQuery(cx,2);
00583   if(fCalibrationSgateToNDNu.GetNumRows()>0) {
00584     static int sCall = 0;
00585     if ( sCall == 0 ) {
00586       MSG("SpillTimeFinder",Msg::kInfo)
00587         << "Initial Spill Timing Calibration constants (SgateToNDNu):" << endl;
00588     }
00589     for (unsigned int i = 0; i<fCalibrationSgateToNDNu.GetNumRows(); i++) {
00590       if ( sCall == 0 ) {
00591         MSG("SpillTimeFinder",Msg::kInfo) 
00592           << Form("%10.5e",fCalibrationSgateToNDNu.GetRow(i)->GetOffset())
00593           << " +- "
00594           << fCalibrationSgateToNDNu.GetRow(i)->GetError() 
00595           << "\t"
00596           << fCalibrationSgateToNDNu.GetRow(i)->GetDescription() << endl;
00597       }
00598       offset += fCalibrationSgateToNDNu.GetRow(i)->GetOffset();
00599     }
00600     sCall++;
00601     return offset;
00602   }
00603 
00604 
00605   // This changes due to changing the fine SGATE delay at the ND.
00606   // Tuning this doesn' affect the offset above, but only changes
00607   // when the neutrinos are seen relative to the SGATE window.
00608    UInt_t sec = cx.GetTimeStamp().GetSec();
00609   if(sec < 1115450000) {
00610     return 0.4 * Munits::microsecond; // 400ns into the window, as reported by Peter Shannahan
00611   } else if(sec < 1115977300) {
00612     return -0.7 * Munits::microsecond; // A mistuning introduced May 6 by me
00613   } else {
00614     return 1.5 * Munits::microsecond; // The final tuning introduced May 12 to give a 1.5us window before the spill.
00615   }
00616 
00617   
00618 }

const SpillTimeND & SpillTimeFinder::GetRecentSpill ( const VldContext context  ) 

Definition at line 686 of file SpillTimeFinder.cxx.

References FindClosestEntries(), and kSpillTime_VeryOld.

Referenced by GetTimeOfRecentSpill().

00687 {
00688   // Return reference to most recent spill record.
00689 
00690   const SpillTimeND* prev;
00691   const SpillTimeND* next;
00692   FindClosestEntries(context,prev,next);
00693   if(prev) return *prev;
00694   else     return kSpillTime_VeryOld;
00695 }

Double_t SpillTimeFinder::GetTimeDifference ( const VldTimeStamp t1,
const VldTimeStamp v2 
) [private]

Definition at line 672 of file SpillTimeFinder.cxx.

References MuELoss::e, VldTimeStamp::GetNanoSec(), and VldTimeStamp::GetSec().

Referenced by FindBestRows(), GetNearestSpill(), GetTimeToNearestSpill(), GetTimeToNextSpill(), and GetTimeToRecentSpill().

00674 {
00675   // Return a-b
00676   double d1 = a.GetSec() - b.GetSec();
00677   double d2 = (a.GetNanoSec() - b.GetNanoSec())*1e-9;
00678   return d1+d2;  
00679 }

VldTimeStamp SpillTimeFinder::GetTimeOfNearestSpill ( const VldContext cx  ) 
VldTimeStamp SpillTimeFinder::GetTimeOfNextSpill ( const VldContext cx  ) 
VldTimeStamp SpillTimeFinder::GetTimeOfRecentSpill ( const VldContext cx  ) 

Definition at line 737 of file SpillTimeFinder.cxx.

References VldTimeStamp::Add(), fTask, GetOffset(), GetRecentSpill(), and SpillTimeND::GetTimeStamp().

Referenced by GetTimeToRecentSpill(), and EventInfoPage::WriteInfo().

00738 {
00739   VldTimeStamp tspill = GetRecentSpill(cx).GetTimeStamp();
00740   tspill.Add(GetOffset(cx,fTask));
00741   return tspill;
00742 }

Double_t SpillTimeFinder::GetTimeToNearestSpill ( const VldContext cx  )  [inline]
Double_t SpillTimeFinder::GetTimeToNextSpill ( const VldContext cx  )  [inline]

Definition at line 110 of file SpillTimeFinder.h.

References GetTimeDifference(), GetTimeOfNextSpill(), and VldContext::GetTimeStamp().

00111 {
00112   return GetTimeDifference(GetTimeOfNextSpill(cx),cx.GetTimeStamp()); 
00113 }

Double_t SpillTimeFinder::GetTimeToRecentSpill ( const VldContext cx  )  [inline]

Definition at line 105 of file SpillTimeFinder.h.

References GetTimeDifference(), GetTimeOfRecentSpill(), and VldContext::GetTimeStamp().

00106 {
00107   return GetTimeDifference(GetTimeOfRecentSpill(cx),cx.GetTimeStamp()); 
00108 }

SpillTimeFinder & SpillTimeFinder::Instance ( void   )  [static]

Instance() returns the singleton handle. ///

Definition at line 15 of file SpillTimeFinder.cxx.

References fgInstance, Msg::kError, MSG, SpillTimeFinder(), and SpillTimeFinder::Cleaner::UseMe().

Referenced by SelectSpillTimes::Ana(), ParticleBeamMonAna::ana(), NDSgateTimeLooter::Ana(), FardetBeamSelect::Ana(), ClockCalibrationModule::Ana(), CheckND::Ana(), BeamMonAna::Analyze(), BDCheckDB::CheckSpill(), CloseToSpillAtFar(), MadTVAnalysis::CreatePAN(), MadTestAnalysis::CreatePAN(), MadPIDAnalysis::CreatePAN(), MadMKAnalysis::CreatePAN(), MadDpAnalysis::CreatePAN(), MadAnalysis::CreatePAN(), BDataQualityModule::EndFile(), NuExtraction::ExtractTimeToNearestSpill(), fill_bdtest(), ANtpInfoObjectFillerBeam::FillBeamInformation(), BDataQualityModule::FillFile(), NtpBDLiteModule::FillNtpBDLite(), NtpMaker::FillSpillInfo(), DataQualityInterface::ProcessBeamStatus(), DataUtil::QueryBeamDB(), NueBeamMonModule::Reco(), Anp::FillHeader::Run(), StndBmsSpin::Scan(), STND_BMS::Scan(), BMS_STND::Scan(), MadQuantities::ShowerValidation(), and EventInfoPage::WriteInfo().

00016 {
00019 
00020   // Cleaner destructor calls DigitSpillTimeFinder dtor
00021   static Cleaner cleaner;
00022   
00023   // SpillTimeFinder mode "1" the default.  This is hardwired for now.
00024   if (!fgInstance) {
00025     cleaner.UseMe();
00026     fgInstance = new SpillTimeFinder();
00027     if (!fgInstance) {
00028       MSG("SpillTime", Msg::kError)
00029         << "No DigitSpillTimeFinder Instance - fatal." << endl;
00030       assert(fgInstance);               // Kill job is there is no instance
00031     }
00032   }
00033   return *fgInstance;
00034 }

void SpillTimeFinder::ResetCache (  )  [private]
static void SpillTimeFinder::SwitchOffsetsOff (  )  [inline, static]

Definition at line 39 of file SpillTimeFinder.h.

References fSwitchOffsetsOff.

00039 {fSwitchOffsetsOff = true;};

static void SpillTimeFinder::SwitchOffsetsOn (  )  [inline, static]

Definition at line 40 of file SpillTimeFinder.h.

References fSwitchOffsetsOff.

00040 {fSwitchOffsetsOff = false;};


Friends And Related Function Documentation

friend struct Cleaner [friend]

Definition at line 97 of file SpillTimeFinder.h.


Member Data Documentation

Definition at line 74 of file SpillTimeFinder.h.

Referenced by GetOffsetKickerToNDNu().

Definition at line 73 of file SpillTimeFinder.h.

Referenced by GetOffsetNDNuToFDNu().

Definition at line 75 of file SpillTimeFinder.h.

Referenced by GetOffsetSgateToNDNu().

Definition at line 69 of file SpillTimeFinder.h.

Referenced by DataIsAvailable(), and FindClosestEntries().

Definition at line 98 of file SpillTimeFinder.h.

Referenced by Instance(), and SpillTimeFinder::Cleaner::~Cleaner().

Definition at line 79 of file SpillTimeFinder.h.

Referenced by FindClosestEntries(), and ResetCache().

Definition at line 78 of file SpillTimeFinder.h.

Referenced by FindClosestEntries(), and ResetCache().

Definition at line 80 of file SpillTimeFinder.h.

Referenced by FindClosestEntries(), and GetLastQueryRangeBeg().

Definition at line 81 of file SpillTimeFinder.h.

Referenced by FindClosestEntries(), and GetLastQueryRangeEnd().

Definition at line 71 of file SpillTimeFinder.h.

Referenced by FindClosestEntries(), and ~SpillTimeFinder().

Definition at line 70 of file SpillTimeFinder.h.

Referenced by FindClosestEntries(), and ~SpillTimeFinder().

Bool_t SpillTimeFinder::fSwitchOffsetsOff = false [static, private]

Definition at line 85 of file SpillTimeFinder.h.

Referenced by GetOffset(), SwitchOffsetsOff(), and SwitchOffsetsOn().

Int_t SpillTimeFinder::fTask [private]

The documentation for this class was generated from the following files:

Generated on 2 Nov 2017 for loon by  doxygen 1.6.1