Data handler computing correlations and linear-regression coefficients. More...

#include <QwCorrelator.h>

Inheritance diagram for QwCorrelator:

Collaboration diagram for QwCorrelator:

Public Member Functions
	QwCorrelator (const TString &name)
	Constructor with name.

	QwCorrelator (const QwCorrelator &name)

	~QwCorrelator () override

void	ProcessOptions (QwOptions &options)

void	ParseConfigFile (QwParameterFile &file) override

Int_t	LoadChannelMap (const std::string &mapfile) override

Int_t	ConnectChannels (QwSubsystemArrayParity &asym, QwSubsystemArrayParity &diff) override
	Connect to Channels (asymmetry/difference only)

void	ProcessData () override

void	FinishDataHandler () override

void	CalcCorrelations ()

void	ConstructTreeBranches (QwRootFile *treerootfile, const std::string &treeprefix="", const std::string &branchprefix="") override
	Construct the tree branches.

void	FillTreeBranches (QwRootFile *treerootfile) override
	Fill the tree branches.

void	ConstructHistograms (TDirectory *folder, TString &prefix) override
	Construct the histograms in a folder with a prefix.

void	FillHistograms () override
	Fill the histograms.

void	ClearEventData () override

void	AccumulateRunningSum (VQwDataHandler &value, Int_t count=0, Int_t ErrorMask=0xFFFFFFF) override

Public Member Functions inherited from VQwDataHandler
	VQwDataHandler (const TString &name)

	VQwDataHandler (const VQwDataHandler &source)

void	SetPointer (QwHelicityPattern *ptr)

void	SetPointer (QwSubsystemArrayParity *ptr)

virtual Int_t	ConnectChannels (QwSubsystemArrayParity &, QwSubsystemArrayParity &asym, QwSubsystemArrayParity &diff)

virtual Int_t	ConnectChannels (QwSubsystemArrayParity &)

Int_t	ConnectChannels (QwHelicityPattern &helicitypattern)

virtual void	UpdateBurstCounter (Short_t burstcounter)

	~VQwDataHandler () override

TString	GetName ()

void	InitRunningSum ()

void	AccumulateRunningSum ()

void	CalculateRunningAverage ()

void	PrintValue () const

void	WritePromptSummary (QwPromptSummary *ps, TString type)

virtual void	ConstructNTupleFields (QwRootFile *treerootfile, const std::string &treeprefix="", const std::string &branchprefix="")
	RNTuple methods.

virtual void	FillNTupleFields (QwRootFile *treerootfile)

void	FillTreeVector (QwRootTreeBranchVector &values) const

void	ConstructBranchAndVector (TTree *tree, TString &prefix, QwRootTreeBranchVector &values)

void	SetRunLabel (TString x)

Int_t	LoadChannelMap ()

Bool_t	PublishInternalValues () const override
	Publish all variables of the subsystem.

Bool_t	PublishByRequest (TString device_name) override
	Try to publish an internal variable matching the submitted name.

Public Member Functions inherited from MQwPublishable_child< QwDataHandlerArray, VQwDataHandler >
	MQwPublishable_child ()
	Default constructor Initializes the child object and sets up self-reference for publishing.

	MQwPublishable_child (const MQwPublishable_child &source)
	Copy constructor.

virtual	~MQwPublishable_child ()
	Virtual destructor.

void	SetParent (QwDataHandlerArray *parent)
	Set the parent container for this child object.

QwDataHandlerArray *	GetParent () const
	Get the parent container for this child object.

Public Member Functions inherited from MQwCloneable< VQwDataHandler, QwCorrelator >
	~MQwCloneable () override
	Virtual destructor.

VQwDataHandler *	Clone () const override
	Concrete clone method.

const VQwFactory< VQwDataHandler > *	Factory () const override
	Factory getter.

Public Member Functions inherited from VQwCloneable< VQwDataHandler >
virtual	~VQwCloneable ()
	Virtual destructor.

std::string	GetClassName () const
	Get demangled name of this class.

Static Public Member Functions
static void	DefineOptions (QwOptions &options)

Static Public Member Functions inherited from MQwCloneable< VQwDataHandler, QwCorrelator >
static VQwDataHandler *	Create (const std::string &name)
	Object creation.

static QwCorrelator *	Cast (QwCorrelator *type)
	Object dynamic cast.

Protected Member Functions
void	OpenAlphaFile (const std::string &prefix)

void	WriteAlphaFile ()

void	CloseAlphaFile ()

void	OpenAliasFile (const std::string &prefix)

void	WriteAliasFile ()

void	CloseAliasFile ()

Protected Member Functions inherited from VQwDataHandler
	VQwDataHandler ()

void	SetEventcutErrorFlagPointer (const UInt_t *errorflagptr)

UInt_t	GetEventcutErrorFlag () const

std::pair< EQwHandleType, std::string >	ParseHandledVariable (const std::string &variable)

void	CalcOneOutput (const VQwHardwareChannel dv, VQwHardwareChannel output, std::vector< const VQwHardwareChannel * > &ivs, std::vector< Double_t > &sens)

Protected Member Functions inherited from MQwPublishable_child< QwDataHandlerArray, VQwDataHandler >
Bool_t	RequestExternalValue (const TString &name, VQwHardwareChannel *value) const
	Retrieve the variable name from other subsystem arrays Get the value corresponding to some variable name from a different data array.

const VQwHardwareChannel *	RequestExternalPointer (const TString &name) const
	Retrieve a pointer to an external variable by name Requests a direct pointer to a variable from sibling subsystems via the parent container.

Bool_t	PublishInternalValue (const TString name, const TString desc, const VQwHardwareChannel *element) const
	Publish a variable from this child into the parent container.

Protected Attributes
Int_t	fBlock

bool	fDisableHistos

std::vector< std::string >	fIndependentFull

std::vector< EQwHandleType >	fIndependentType

std::vector< std::string >	fIndependentName

std::vector< const VQwHardwareChannel * >	fIndependentVar

std::vector< Double_t >	fIndependentValues

std::string	fAlphaOutputFileBase

std::string	fAlphaOutputFileSuff

std::string	fAlphaOutputPath

TFile *	fAlphaOutputFile

TTree *	fTree

std::string	fAliasOutputFileBase

std::string	fAliasOutputFileSuff

std::string	fAliasOutputPath

std::ofstream	fAliasOutputFile

int	fTotalCount

int	fGoodCount

int	fErrCounts_EF

std::vector< int >	fErrCounts_IV

std::vector< int >	fErrCounts_DV

unsigned int	fGoodEvent

Protected Attributes inherited from VQwDataHandler
Int_t	fPriority

Short_t	fBurstCounter
	When a datahandler array is processed, handlers with lower priority will be processed before handlers with higher priority.

TString	fName

std::string	fMapFile

std::string	fTreeName

std::string	fTreeComment

std::string	fPrefix

TString	run_label

const UInt_t *	fErrorFlagPtr
	Error flag pointer.

QwSubsystemArrayParity *	fSubsystemArray
	Single event pointer.

QwHelicityPattern *	fHelicityPattern
	Helicity pattern pointer.

std::vector< std::string >	fDependentFull

std::vector< EQwHandleType >	fDependentType

std::vector< std::string >	fDependentName

std::vector< const VQwHardwareChannel * >	fDependentVar

std::vector< Double_t >	fDependentValues

std::vector< VQwHardwareChannel * >	fOutputVar

std::vector< Double_t >	fOutputValues

std::vector< std::vector< TString > >	fPublishList

std::string	ParseSeparator

Bool_t	fKeepRunningSum

Bool_t	fRunningsumFillsTree

VQwDataHandler *	fRunningsum

Private Member Functions
	QwCorrelator ()

Private Attributes
TString	fNameNoSpaces

int	nP

int	nY

std::vector< TH1D >	fHnames

std::vector< TH1D >	fH1iv

std::vector< TH1D >	fH1dv

std::vector< std::vector< TH2D > >	fH2iv

std::vector< std::vector< TH2D > >	fH2dv

LinRegBevPeb	linReg

Int_t	fCycleCounter

Static Private Attributes
static bool	fPrintCorrelations = false

Additional Inherited Members
Public Types inherited from VQwDataHandler
enum	EQwHandleType { kHandleTypeUnknown =0 , kHandleTypeMps , kHandleTypeAsym , kHandleTypeDiff , kHandleTypeYield }

typedef std::vector< VQwHardwareChannel * >::iterator	Iterator_HdwChan

typedef std::vector< VQwHardwareChannel * >::const_iterator	ConstIterator_HdwChan

Detailed Description

Data handler computing correlations and linear-regression coefficients.

Uses Bevington/Pebay algorithms to estimate correlations between independent and dependent variables selected from subsystem arrays. Produces summary histograms and optional output trees/files for further analysis.

Definition at line 33 of file QwCorrelator.h.

Constructor & Destructor Documentation

◆ QwCorrelator() [1/3]

QwCorrelator::QwCorrelator ( const TString & name )

Constructor with name.

Definition at line 28 of file QwCorrelator.cc.

: VQwDataHandler(name),
  fBlock(-1),
  fDisableHistos(true),
  fAlphaOutputFileBase("blueR"),
  fAlphaOutputFileSuff("new.slope.root"),
  fAlphaOutputPath("."),
  fAlphaOutputFile(0),
  fTree(0),
  fAliasOutputFileBase("regalias_"),
  fAliasOutputFileSuff(""),
  fAliasOutputPath("."),
  fNameNoSpaces(name),
  nP(0),nY(0),
  fCycleCounter(0)
{
  fNameNoSpaces.ReplaceAll(" ","_");
  // Set default tree name and descriptions (in VQwDataHandler)
  fTreeName = "lrb";
  fTreeComment = "Correlations";
  // Parsing separator
  ParseSeparator = "_";
 
  // Clear all data
  ClearEventData();
}

References ClearEventData(), fAliasOutputFileBase, fAliasOutputFileSuff, fAliasOutputPath, fAlphaOutputFile, fAlphaOutputFileBase, fAlphaOutputFileSuff, fAlphaOutputPath, fBlock, fCycleCounter, fDisableHistos, fNameNoSpaces, fTree, VQwDataHandler::fTreeComment, VQwDataHandler::fTreeName, nP, nY, VQwDataHandler::ParseSeparator, and VQwDataHandler::VQwDataHandler().

Referenced by AccumulateRunningSum(), and QwCorrelator().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ QwCorrelator() [2/3]

QwCorrelator::QwCorrelator ( const QwCorrelator & name )

Definition at line 55 of file QwCorrelator.cc.

: VQwDataHandler(source),
  fBlock(source.fBlock),
  fDisableHistos(source.fDisableHistos),
  fAlphaOutputFileBase(source.fAlphaOutputFileBase),
  fAlphaOutputFileSuff(source.fAlphaOutputFileSuff),
  fAlphaOutputPath(source.fAlphaOutputPath),
  fAlphaOutputFile(0),
  fTree(0),
  fAliasOutputFileBase(source.fAliasOutputFileBase),
  fAliasOutputFileSuff(source.fAliasOutputFileSuff),
  fAliasOutputPath(source.fAliasOutputPath),
  nP(source.nP),nY(source.nY),
  fCycleCounter(source.fCycleCounter)
{
  QwWarning << "QwCorrelator copy constructor required but untested" << QwLog::endl;
 
  // Clear all data
  ClearEventData();
}

References ClearEventData(), QwLog::endl(), fAliasOutputFileBase, fAliasOutputFileSuff, fAliasOutputPath, fAlphaOutputFile, fAlphaOutputFileBase, fAlphaOutputFileSuff, fAlphaOutputPath, fBlock, fCycleCounter, fDisableHistos, fTree, nP, nY, QwCorrelator(), QwWarning, and VQwDataHandler::VQwDataHandler().

Here is the call graph for this function:

◆ ~QwCorrelator()

QwCorrelator::~QwCorrelator ( )

override

Definition at line 76 of file QwCorrelator.cc.

{
  // Close alpha and alias file
  CloseAlphaFile();
  CloseAliasFile();
}

References CloseAliasFile(), and CloseAlphaFile().

Here is the call graph for this function:

◆ QwCorrelator() [3/3]

QwCorrelator::QwCorrelator ( )

private

Member Function Documentation

◆ AccumulateRunningSum()

void QwCorrelator::AccumulateRunningSum	(	VQwDataHandler &	value,
		Int_t	count = 0,
		Int_t	ErrorMask = 0xFFFFFFF )

overridevirtual

Reimplemented from VQwDataHandler.

Definition at line 160 of file QwCorrelator.cc.

{
  QwCorrelator* correlator = dynamic_cast<QwCorrelator*>(&value);
  if (correlator) {
    linReg += correlator->linReg;
  } else {
    QwWarning << "QwCorrelator::AccumulateRunningSum "
              << "can only accept other QwCorrelator objects."
              << QwLog::endl;
  }
}

References QwLog::endl(), linReg, QwCorrelator(), QwWarning, and VQwDataHandler::VQwDataHandler().

Here is the call graph for this function:

◆ CalcCorrelations()

void QwCorrelator::CalcCorrelations ( )

Definition at line 172 of file QwCorrelator.cc.

{
  // Check if any channels are active
  if (nP == 0 || nY == 0) {
    return;
  }
 
  QwMessage << "QwCorrelator::CalcCorrelations(): name=" << GetName() << QwLog::endl;
 
  // Print entry summary
  QwVerbose << "QwCorrelator: "
            << "total entries: " << fTotalCount << ", "
            << "good entries: " << fGoodCount
            << QwLog::endl;
  // and warn if zero
  if (fTotalCount > 100 && fGoodCount == 0) {
    QwWarning << "QwCorrelator: "
              << "< 1% good events, "
              << fGoodCount << " of " << fTotalCount
              << QwLog::endl;
  }
 
  // Event error flag
  if (fErrCounts_EF > 0) {
    QwVerbose << "   Entries failed due to error flag: "
              << fErrCounts_EF << QwLog::endl;
  }
  // Dependent variable error codes
  for (size_t i = 0; i < fDependentVar.size(); ++i) {
    if (fErrCounts_DV.at(i) > 0) {
      QwVerbose << "   Entries failed due to " << fDependentVar.at(i)->GetElementName()
                << ": " <<  fErrCounts_DV.at(i) << QwLog::endl;
    }
  }
  // Independent variable error codes
  for (size_t i = 0; i < fIndependentVar.size(); ++i) {
    if (fErrCounts_IV.at(i) > 0) {
      QwVerbose << "   Entries failed due to " << fIndependentVar.at(i)->GetElementName()
                << ": " <<  fErrCounts_IV.at(i) << QwLog::endl;
    }
  }
 
  if (! linReg.failed()) {
 
    if (fPrintCorrelations) {
      linReg.printSummaryP();
      linReg.printSummaryY();
    }
 
    linReg.solve();
 
    if (fPrintCorrelations) {
      linReg.printSummaryAlphas();
      linReg.printSummaryMeansWithUnc();
      linReg.printSummaryMeansWithUncCorrected();
    }
  }
 
  // Fill tree
  if (fTree) fTree->Fill();
  else QwWarning << "No tree" << QwLog::endl;
 
  // Write alpha and alias file
  WriteAlphaFile();
  WriteAliasFile();
}

References QwLog::endl(), VQwDataHandler::fDependentVar, fErrCounts_DV, fErrCounts_EF, fErrCounts_IV, fGoodCount, fIndependentVar, fPrintCorrelations, fTotalCount, fTree, VQwDataHandler::GetName(), linReg, nP, nY, QwMessage, QwVerbose, QwWarning, WriteAliasFile(), and WriteAlphaFile().

Referenced by FinishDataHandler().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ ClearEventData()

void QwCorrelator::ClearEventData ( )

overridevirtual

Reimplemented from VQwDataHandler.

Definition at line 144 of file QwCorrelator.cc.

{
  // Clear error counters
  fErrCounts_EF = 0;
  std::fill(fErrCounts_DV.begin(), fErrCounts_DV.end(), 0);
  std::fill(fErrCounts_IV.begin(), fErrCounts_IV.end(), 0);
 
  // Clear event counts
  fTotalCount = 0;
  fGoodCount = 0;
  fGoodEvent = -1;
 
  // Clear regression
  linReg.clear();
}

References fErrCounts_DV, fErrCounts_EF, fErrCounts_IV, fGoodCount, fGoodEvent, fTotalCount, and linReg.

Referenced by QwCorrelator(), and QwCorrelator().

Here is the caller graph for this function:

◆ CloseAliasFile()

void QwCorrelator::CloseAliasFile ( )

protected

Definition at line 688 of file QwCorrelator.cc.

{
  // Close alias output file
  if (fAliasOutputFile.good()) {
    fAliasOutputFile << "}" << std::endl << std::endl;
    fAliasOutputFile.close();
  } else {
    QwWarning << "QwCorrelator: Unable to close alias output file." << QwLog::endl;
  }
}

References QwLog::endl(), fAliasOutputFile, and QwWarning.

Referenced by ~QwCorrelator().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ CloseAlphaFile()

void QwCorrelator::CloseAlphaFile ( )

protected

Definition at line 676 of file QwCorrelator.cc.

{
  // Close slopes output file
  if (fAlphaOutputFile) {
    // Switch to this file's directory before writing to avoid
    // trying to write objects from other (possibly closed) files
    fAlphaOutputFile->cd();
    fAlphaOutputFile->Write();
    fAlphaOutputFile->Close();
  }
}

References fAlphaOutputFile.

Referenced by ~QwCorrelator().

Here is the caller graph for this function:

◆ ConnectChannels()

Int_t QwCorrelator::ConnectChannels	(	QwSubsystemArrayParity &	asym,
		QwSubsystemArrayParity &	diff )

overridevirtual

Connect to Channels (asymmetry/difference only)

Parameters

asym	Subsystem array with asymmetries
diff	Subsystem array with differences

Returns: 0 on success, non-zero on failure

Fill vector of pointers to the relevant data elements

Reimplemented from VQwDataHandler.

Definition at line 287 of file QwCorrelator.cc.

{
  SetEventcutErrorFlagPointer(asym.GetEventcutErrorFlagPointer());
 
  // Return if correlator is not enabled
 
  /// Fill vector of pointers to the relevant data elements
  for (size_t dv = 0; dv < fDependentName.size(); dv++) {
    // Get the dependent variables
 
    const VQwHardwareChannel* dv_ptr = 0;
 
    if (fDependentType.at(dv)==kHandleTypeMps){
      //  Quietly ignore the MPS type when we're connecting the asym & diff
      continue;
    }else{
      dv_ptr = this->RequestExternalPointer(fDependentFull.at(dv));
      if (dv_ptr==NULL){
    switch (fDependentType.at(dv)) {
        case kHandleTypeAsym:
          dv_ptr = asym.RequestExternalPointer(fDependentName.at(dv));
          break;
        case kHandleTypeDiff:
          dv_ptr = diff.RequestExternalPointer(fDependentName.at(dv));
          break;
        default:
          QwWarning << "QwCorrelator::ConnectChannels(QwSubsystemArrayParity& asym, QwSubsystemArrayParity& diff): "
                    << "Dependent variable, " << fDependentName.at(dv)
                    << ", for asym/diff correlator does not have proper type, type=="
                    << fDependentType.at(dv) << "." << QwLog::endl;
          break;
        }
      }
      if (dv_ptr == NULL){
    QwWarning << "QwCombiner::ConnectChannels(QwSubsystemArrayParity& asym, QwSubsystemArrayParity& diff):  Dependent variable, "
          << fDependentName.at(dv)
          << ", was not found (fullname=="
          << fDependentFull.at(dv)<< ")." << QwLog::endl;
     continue;
      }
    }
 
    // pair creation
    if(dv_ptr != NULL){
      // fDependentVarType.push_back(fDependentType.at(dv));
      fDependentVar.push_back(dv_ptr);
    }
 
  }
 
  // Add independent variables
  for (size_t iv = 0; iv < fIndependentName.size(); iv++) {
    // Get the independent variables
    const VQwHardwareChannel* iv_ptr = 0;
    iv_ptr = this->RequestExternalPointer(fIndependentFull.at(iv));
    if (iv_ptr==NULL){
      switch (fIndependentType.at(iv)) {
      case kHandleTypeAsym:
        iv_ptr = asym.RequestExternalPointer(fIndependentName.at(iv));
        break;
      case kHandleTypeDiff:
        iv_ptr = diff.RequestExternalPointer(fIndependentName.at(iv));
        break;
      default:
        QwWarning << "Independent variable for correlator has unknown type."
                  << QwLog::endl;
        break;
      }
    }
    if (iv_ptr) {
      fIndependentVar.push_back(iv_ptr);
    } else {
      QwWarning << "Independent variable " << fIndependentName.at(iv) << " for correlator could not be found."
                << QwLog::endl;
    }
 
  }
  fIndependentValues.resize(fIndependentVar.size());
  fDependentValues.resize(fDependentVar.size());
 
  nP = fIndependentName.size();
  nY = fDependentName.size();
 
  linReg.setDims(nP, nY);
  linReg.init();
 
  fErrCounts_IV.resize(fIndependentVar.size(),0);
  fErrCounts_DV.resize(fDependentVar.size(),0);
 
  return 0;
}

References QwLog::endl(), VQwDataHandler::fDependentFull, VQwDataHandler::fDependentName, VQwDataHandler::fDependentType, VQwDataHandler::fDependentValues, VQwDataHandler::fDependentVar, fErrCounts_DV, fErrCounts_IV, fIndependentFull, fIndependentName, fIndependentType, fIndependentValues, fIndependentVar, QwSubsystemArrayParity::GetEventcutErrorFlagPointer(), VQwDataHandler::kHandleTypeAsym, VQwDataHandler::kHandleTypeDiff, VQwDataHandler::kHandleTypeMps, linReg, nP, nY, QwWarning, MQwPublishable< U, T >::RequestExternalPointer(), MQwPublishable_child< QwDataHandlerArray, VQwDataHandler >::RequestExternalPointer(), and VQwDataHandler::SetEventcutErrorFlagPointer().

Here is the call graph for this function:

◆ ConstructHistograms()

void QwCorrelator::ConstructHistograms	(	TDirectory *	folder,
		TString &	prefix )

overridevirtual

Construct the histograms in a folder with a prefix.

Reimplemented from VQwDataHandler.

Definition at line 468 of file QwCorrelator.cc.

{
  // Skip if disabled
  if (fDisableHistos) return;
 
  // Check if any channels are active
  if (nP == 0 || nY == 0) {
    return;
  }
 
  // Go to directory
  TString name(fName);
  name.ReplaceAll(" ","_");
  folder->mkdir(name)->cd();
 
  //..... 1D,  iv
  fH1iv.resize(nP);
  for (int i = 0; i < nP; i++) {
    fH1iv[i] = TH1D(
        Form("P%d",i),
        Form("iv P%d=%s, pass=%s ;iv=%s (ppm)",i,fIndependentName[i].c_str(),fName.Data(),fIndependentName[i].c_str()),
        128,0.,0.);
    fH1iv[i].GetXaxis()->SetNdivisions(4);
  }
 
  //..... 2D,  iv correlations
  Double_t x1 = 0;
  fH2iv.resize(nP);
  for (int i = 0; i < nP; i++) {
    fH2iv[i].resize(nP);
    for (int j = i+1; j < nP; j++) { // not all are used
      fH2iv[i][j] = TH2D(
          Form("P%d_P%d",i,j),
          Form("iv correlation  P%d_P%d, pass=%s ;P%d=%s (ppm);P%d=%s   (ppm)  ",
              i,j,fName.Data(),i,fIndependentName[i].c_str(),j,fIndependentName[j].c_str()),
          64,-x1,x1,
          64,-x1,x1);
      fH2iv[i][j].GetXaxis()->SetTitleColor(kBlue);
      fH2iv[i][j].GetYaxis()->SetTitleColor(kBlue);
      fH2iv[i][j].GetXaxis()->SetNdivisions(4);
      fH2iv[i][j].GetYaxis()->SetNdivisions(4);
    }
  }
 
  //..... 1D,  dv
  fH1dv.resize(nY);
  for (int i = 0; i < nY; i++) {
    fH1dv[i] = TH1D(
        Form("Y%d",i),
        Form("dv Y%d=%s, pass=%s ;dv=%s (ppm)",i,fDependentName[i].c_str(),fName.Data(),fDependentName[i].c_str()),
        128,0.,0.);
    fH1dv[i].GetXaxis()->SetNdivisions(4);
  }
 
  //..... 2D,  dv-iv correlations
  Double_t y1 = 0;
  fH2dv.resize(nP);
  for (int i = 0; i < nP; i++) {
    fH2dv[i].resize(nY);
    for (int j = 0; j < nY; j++) {
      fH2dv[i][j] = TH2D(
          Form("P%d_Y%d",i,j),
          Form("iv-dv correlation  P%d_Y%d, pass=%s ;P%d=%s (ppm);Y%d=%s   (ppm)  ",
              i,j,fName.Data(),i,fIndependentName[i].c_str(),j,fDependentName[j].c_str()),
          64,-x1,x1,
          64,-y1,y1);
      fH2dv[i][j].GetXaxis()->SetTitleColor(kBlue);
      fH2dv[i][j].GetYaxis()->SetTitleColor(kBlue);
      fH2dv[i][j].GetXaxis()->SetNdivisions(4);
      fH2dv[i][j].GetYaxis()->SetNdivisions(4);
    }
  }
 
  // store list of names to be archived
  fHnames.resize(2);
  fHnames[0] = TH1D("NamesIV",Form("IV name list nIV=%d",nP),nP,0,1);
  for (int i = 0; i < nP; i++)
    fHnames[0].Fill(fIndependentName[i].c_str(),1.*i);
  fHnames[1] = TH1D("NamesDV",Form("DV name list nIV=%d",nY),nY,0,1);
  for (int i = 0; i < nY; i++)
    fHnames[1].Fill(fDependentName[i].c_str(),i*1.);
}

References VQwDataHandler::fDependentName, fDisableHistos, fH1dv, fH1iv, fH2dv, fH2iv, fHnames, fIndependentName, VQwDataHandler::fName, nP, and nY.

◆ ConstructTreeBranches()

void QwCorrelator::ConstructTreeBranches	(	QwRootFile *	treerootfile,
		const std::string &	treeprefix = "",
		const std::string &	branchprefix = "" )

overridevirtual

Construct the tree branches.

Reimplemented from VQwDataHandler.

Definition at line 380 of file QwCorrelator.cc.

{
  // Check if any channels are active
  if (nP == 0 || nY == 0) {
    return;
  }
 
  // Check if tree name is specified
  if (fTreeName == "") {
    QwWarning << "QwCorrelator: no tree name specified, use 'tree-name = value'" << QwLog::endl;
    return;
  }
 
  // Create alpha and alias files before trying to create the tree
  OpenAlphaFile(treeprefix);
  OpenAliasFile(treeprefix);
 
  // Construct tree name and create new tree
  const std::string name = treeprefix + fTreeName;
  treerootfile->NewTree(name, fTreeComment.c_str());
  fTree = treerootfile->GetTree(name);
  // Check to make sure the tree was created successfully
  if (fTree == NULL) return;
 
  // Set up branches
  fTree->Branch(TString(branchprefix + "total_count"), &fTotalCount);
  fTree->Branch(TString(branchprefix + "good_count"),  &fGoodCount);
 
  fTree->Branch(TString(branchprefix + "n"), &(linReg.fGoodEventNumber));
  fTree->Branch(TString(branchprefix + "ErrorFlag"), &(linReg.fErrorFlag));
 
  auto bn = [&](const TString& n) {
    return TString(branchprefix + n);
  };
  auto pm = [](TMatrixD& m) {
    return m.GetMatrixArray();
  };
  auto lm = [](TMatrixD& m, const TString& n) {
    return Form("%s[%d][%d]/D", n.Data(), m.GetNrows(), m.GetNcols());
  };
  auto branchm = [&](TTree* tree, TMatrixD& m, const TString& n) {
    tree->Branch(bn(n),pm(m),lm(m,n));
  };
  auto pv = [](TVectorD& v) {
    return v.GetMatrixArray();
  };
  auto lv = [](TVectorD& v, const TString& n) {
    return Form("%s[%d]/D", n.Data(), v.GetNrows());
  };
  auto branchv = [&](TTree* tree, TVectorD& v, const TString& n) {
    tree->Branch(bn(n),pv(v),lv(v,n));
  };
 
  branchm(fTree,linReg.Axy,  "A");
  branchm(fTree,linReg.dAxy, "dA");
 
  branchm(fTree,linReg.mVPP,  "VPP");
  branchm(fTree,linReg.mVPY,  "VPY");
  branchm(fTree,linReg.mVYP,  "VYP");
  branchm(fTree,linReg.mVYY,  "VYY");
  branchm(fTree,linReg.mVYYp, "VYYp");
 
  branchm(fTree,linReg.mSPP,  "SPP");
  branchm(fTree,linReg.mSPY,  "SPY");
  branchm(fTree,linReg.mSYP,  "SYP");
  branchm(fTree,linReg.mSYY,  "SYY");
  branchm(fTree,linReg.mSYYp, "SYYp");
 
  branchm(fTree,linReg.mRPP,  "RPP");
  branchm(fTree,linReg.mRPY,  "RPY");
  branchm(fTree,linReg.mRYP,  "RYP");
  branchm(fTree,linReg.mRYY,  "RYY");
  branchm(fTree,linReg.mRYYp, "RYYp");
 
  branchv(fTree,linReg.mMP,  "MP");   // Parameter mean
  branchv(fTree,linReg.mMY,  "MY");   // Uncorrected mean
  branchv(fTree,linReg.mMYp, "MYp");  // Corrected mean
 
  branchv(fTree,linReg.mSP,  "dMP");  // Parameter mean error
  branchv(fTree,linReg.mSY,  "dMY");  // Uncorrected mean error
  branchv(fTree,linReg.mSYp, "dMYp"); // Corrected mean error
 
}

References QwLog::endl(), fGoodCount, fTotalCount, fTree, VQwDataHandler::fTreeComment, VQwDataHandler::fTreeName, QwRootFile::GetTree(), linReg, QwRootFile::NewTree(), nP, nY, OpenAliasFile(), OpenAlphaFile(), and QwWarning.

Here is the call graph for this function:

◆ DefineOptions()

void QwCorrelator::DefineOptions ( QwOptions & options )

static

Definition at line 83 of file QwCorrelator.cc.

{
  options.AddOptions()("print-correlations",
      po::value<bool>(&fPrintCorrelations)->default_bool_value(false),
      "print correlations after determining them");
}

References QwOptions::AddOptions(), and fPrintCorrelations.

Referenced by DefineOptionsParity().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ FillHistograms()

void QwCorrelator::FillHistograms ( )

overridevirtual

Fill the histograms.

Reimplemented from VQwDataHandler.

Definition at line 552 of file QwCorrelator.cc.

{
  // Skip if disabled
  if (fDisableHistos) return;
 
  // Check if any channels are active
  if (nP == 0 || nY == 0) {
    return;
  }
 
  // Skip if bad event
  if (fGoodEvent != 0) return;
 
  // Fill histograms
  for (size_t i = 0; i < fIndependentValues.size(); i++) {
    fH1iv[i].Fill(fIndependentValues[i]);
    for (size_t j = i+1; j < fIndependentValues.size(); j++)
      fH2iv[i][j].Fill(fIndependentValues[i], fIndependentValues[j]);
  }
  for (size_t j = 0; j < fDependentValues.size(); j++) {
    fH1dv[j].Fill(fDependentValues[j]);
    for (size_t i = 0; i < fIndependentValues.size(); i++)
      fH2dv[i][j].Fill(fIndependentValues[i], fDependentValues[j]);
  }
}

References VQwDataHandler::fDependentValues, fDisableHistos, fGoodEvent, fH1dv, fH1iv, fH2dv, fH2iv, fIndependentValues, nP, and nY.

◆ FillTreeBranches()

void QwCorrelator::FillTreeBranches ( QwRootFile * treerootfile )

inlineoverridevirtual

Fill the tree branches.

Reimplemented from VQwDataHandler.

Definition at line 69 of file QwCorrelator.h.

69{ };

◆ FinishDataHandler()

void QwCorrelator::FinishDataHandler ( )

inlineoverridevirtual

Reimplemented from VQwDataHandler.

Definition at line 58 of file QwCorrelator.h.

                                   {
    CalcCorrelations();
  }

References CalcCorrelations().

Here is the call graph for this function:

◆ LoadChannelMap()

Int_t QwCorrelator::LoadChannelMap ( const std::string & mapfile )

overridevirtual

Load the channel map

Parameters

mapfile Filename of map file

Returns: Zero when success

Reimplemented from VQwDataHandler.

Definition at line 245 of file QwCorrelator.cc.

{
  // Open the file
  QwParameterFile map(mapfile);
 
  // Read the sections of dependent variables
  std::pair<EQwHandleType,std::string> type_name;
 
  // Add independent variables and sensitivities
  while (map.ReadNextLine()) {
    // Throw away comments, whitespace, empty lines
    map.TrimComment();
    map.TrimWhitespace();
    if (map.LineIsEmpty()) continue;
    // Get first token: label (dv or iv), second token is the name like "asym_blah"
    string primary_token = map.GetNextToken(" ");
    string current_token = map.GetNextToken(" ");
    // Parse current token into independent variable type and name
    type_name = ParseHandledVariable(current_token);
 
    if (primary_token == "iv") {
      fIndependentType.push_back(type_name.first);
      fIndependentName.push_back(type_name.second);
      fIndependentFull.push_back(current_token);
    }
    else if (primary_token == "dv") {
      fDependentType.push_back(type_name.first);
      fDependentName.push_back(type_name.second);
      fDependentFull.push_back(current_token);
    }
    else if (primary_token == "treetype") {
      QwMessage << "Tree Type read, ignoring." << QwLog::endl;
    }
    else {
      QwError << "LoadChannelMap in QwCorrelator read invalid primary_token " << primary_token << QwLog::endl;
    }
  }
 
  return 0;
}

References QwLog::endl(), VQwDataHandler::fDependentFull, VQwDataHandler::fDependentName, VQwDataHandler::fDependentType, fIndependentFull, fIndependentName, fIndependentType, QwParameterFile::GetNextToken(), QwParameterFile::LineIsEmpty(), VQwDataHandler::ParseHandledVariable(), QwError, QwMessage, QwParameterFile::ReadNextLine(), QwParameterFile::TrimComment(), and QwParameterFile::TrimWhitespace().

Here is the call graph for this function:

◆ OpenAliasFile()

void QwCorrelator::OpenAliasFile ( const std::string & prefix )

protected

Definition at line 658 of file QwCorrelator.cc.

{
  // Turn "." into "_" in run_label (no "." allowed in function name, and must
  // agree with the filename)
  std::string label(run_label);
  std::replace(label.begin(), label.end(), '.', '_');
  // Create old-style regalias script
  std::string name = prefix + fAliasOutputFileBase + label + fAliasOutputFileSuff;
  std::string path = fAliasOutputPath + "/";
  std::string file = path + name + ".C"; // add extension outside of file suffix
  fAliasOutputFile.open(file, std::ofstream::out);
  if (fAliasOutputFile.good()) {
    fAliasOutputFile << Form("void %s(int i = 0) {", name.c_str()) << std::endl;
  } else {
    QwWarning << "QwCorrelator: Could not write to alias output file " << QwLog::endl;
  }
}

References QwLog::endl(), fAliasOutputFile, fAliasOutputFileBase, fAliasOutputFileSuff, fAliasOutputPath, QwWarning, and VQwDataHandler::run_label.

Referenced by ConstructTreeBranches().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ OpenAlphaFile()

void QwCorrelator::OpenAlphaFile ( const std::string & prefix )

protected

Definition at line 644 of file QwCorrelator.cc.

{
  // Create old-style blueR ROOT file
  std::string name = prefix + fAlphaOutputFileBase + run_label.Data() + fAlphaOutputFileSuff;
  std::string path = fAlphaOutputPath + "/";
  std::string file = path + name;
  fAlphaOutputFile = new TFile(TString(file), "RECREATE", "correlation coefficients");
  if (! fAlphaOutputFile->IsWritable()) {
    QwError << "QwCorrelator could not create output file " << file << QwLog::endl;
    delete fAlphaOutputFile;
    fAlphaOutputFile = 0;
  }
}

References QwLog::endl(), fAlphaOutputFile, fAlphaOutputFileBase, fAlphaOutputFileSuff, fAlphaOutputPath, QwError, and VQwDataHandler::run_label.

Referenced by ConstructTreeBranches().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ ParseConfigFile()

void QwCorrelator::ParseConfigFile ( QwParameterFile & file )

overridevirtual

Parse configuration file for map file, priority, tree settings.

Reimplemented from VQwDataHandler.

Definition at line 94 of file QwCorrelator.cc.

{
  VQwDataHandler::ParseConfigFile(file);
  file.PopValue("slope-file-base", fAlphaOutputFileBase);
  file.PopValue("slope-file-suff", fAlphaOutputFileSuff);
  file.PopValue("slope-path", fAlphaOutputPath);
  file.PopValue("alias-file-base", fAliasOutputFileBase);
  file.PopValue("alias-file-suff", fAliasOutputFileSuff);
  file.PopValue("alias-path", fAliasOutputPath);
  file.PopValue("disable-histos", fDisableHistos);
  file.PopValue("block", fBlock);
  if (fBlock >= 4)
    QwWarning << "QwCorrelator: expect 0 <= block <= 3 but block = "
              << fBlock << QwLog::endl;
}

References QwLog::endl(), fAliasOutputFileBase, fAliasOutputFileSuff, fAliasOutputPath, fAlphaOutputFileBase, fAlphaOutputFileSuff, fAlphaOutputPath, fBlock, fDisableHistos, VQwDataHandler::ParseConfigFile(), QwParameterFile::PopValue(), and QwWarning.

Here is the call graph for this function:

◆ ProcessData()

void QwCorrelator::ProcessData ( )

overridevirtual

Copy dependent variables to output variables (default processing).

Reimplemented from VQwDataHandler.

Definition at line 110 of file QwCorrelator.cc.

{
  // Add to total count
  fTotalCount++;
 
  // Start as good event
  fGoodEvent = 0;
 
  // Event error flag
  fGoodEvent |= GetEventcutErrorFlag();
  if ( GetEventcutErrorFlag() != 0) fErrCounts_EF++;
  // Dependent variable error codes
  for (size_t i = 0; i < fDependentVar.size(); ++i) {
    fGoodEvent |= fDependentVar.at(i)->GetErrorCode();
    fDependentValues.at(i) = (fDependentVar[i]->GetValue(fBlock+1));
    if (fDependentVar.at(i)->GetErrorCode() !=0)  (fErrCounts_DV.at(i))++;
  }
  // Independent variable error codes
  for (size_t i = 0; i < fIndependentVar.size(); ++i) {
    fGoodEvent |= fIndependentVar.at(i)->GetErrorCode();
    fIndependentValues.at(i) = (fIndependentVar[i]->GetValue(fBlock+1));
    if (fIndependentVar.at(i)->GetErrorCode() !=0)  (fErrCounts_IV.at(i))++;
  }
 
  // If good, process event
  if (fGoodEvent == 0) {
    fGoodCount++;
 
    TVectorD P(fIndependentValues.size(), fIndependentValues.data());
    TVectorD Y(fDependentValues.size(),   fDependentValues.data());
    linReg += std::make_pair(P, Y);
  }
}

References fBlock, VQwDataHandler::fDependentValues, VQwDataHandler::fDependentVar, fErrCounts_DV, fErrCounts_EF, fErrCounts_IV, fGoodCount, fGoodEvent, fIndependentValues, fIndependentVar, fTotalCount, VQwDataHandler::GetEventcutErrorFlag(), and linReg.

Here is the call graph for this function:

◆ ProcessOptions()

void QwCorrelator::ProcessOptions ( QwOptions & options )

Definition at line 90 of file QwCorrelator.cc.

91{

92}

◆ WriteAliasFile()

void QwCorrelator::WriteAliasFile ( )

protected

Definition at line 699 of file QwCorrelator.cc.

{
  // Ensure output file is open
  if (fAliasOutputFile.bad()) {
    QwWarning << "QwCorrelator: Could not write to alias output file " << QwLog::endl;
    return;
  }
 
  fAliasOutputFile << " if (i == " << fCycleCounter << ") {" << std::endl;
  fAliasOutputFile << Form("  TTree* tree = (TTree*) gDirectory->Get(\"mul\");") << std::endl;
  for (int i = 0; i < nY; i++) {
    fAliasOutputFile << Form("  tree->SetAlias(\"reg_%s\",",fDependentFull[i].c_str()) << std::endl;
    fAliasOutputFile << Form("         \"%s",fDependentFull[i].c_str());
    for (int j = 0; j < nP; j++) {
      fAliasOutputFile << Form("%+.4e*%s", -linReg.Axy(j,i), fIndependentFull[j].c_str());
    }
    fAliasOutputFile << "\");" << std::endl;
  }
  fAliasOutputFile << " }" << std::endl;
 
  // Increment call counter
  fCycleCounter++;
}

References QwLog::endl(), fAliasOutputFile, fCycleCounter, VQwDataHandler::fDependentFull, fIndependentFull, linReg, nP, nY, and QwWarning.

Referenced by CalcCorrelations().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ WriteAlphaFile()

void QwCorrelator::WriteAlphaFile ( )

protected

Definition at line 578 of file QwCorrelator.cc.

{
  // Ensure in output file
  if (fAlphaOutputFile) fAlphaOutputFile->cd();
 
  // Write objects
  linReg.Axy.Write("slopes");
  linReg.dAxy.Write("sigSlopes");
 
  linReg.mRPP.Write("IV_IV_correlation");
  linReg.mRPY.Write("IV_DV_correlation");
  linReg.mRYY.Write("DV_DV_correlation");
  linReg.mRYYp.Write("DV_DV_correlation_prime");
 
  linReg.mMP.Write("IV_mean");
  linReg.mMY.Write("DV_mean");
  linReg.mMYp.Write("DV_mean_prime");
 
  // number of events
  TMatrixD Mstat(1,1);
  Mstat(0,0)=linReg.getUsedEve();
  Mstat.Write("MyStat");
 
  //... IVs
  TH1D hiv("IVname","names of IVs",nP,-0.5,nP-0.5);
  for (int i=0;i<nP;i++) hiv.Fill(fIndependentFull[i].c_str(),i);
  hiv.Write();
 
  //... DVs
  TH1D hdv("DVname","names of IVs",nY,-0.5,nY-0.5);
  for (int i=0;i<nY;i++) hdv.Fill(fDependentFull[i].c_str(),i);
  hdv.Write();
 
  // sigmas
  linReg.mSP.Write("IV_sigma");
  linReg.mSY.Write("DV_sigma");
  linReg.mSYp.Write("DV_sigma_prime");
 
  // raw covariances
  linReg.mVPP.Write("IV_IV_rawVariance");
  linReg.mVPY.Write("IV_DV_rawVariance");
  linReg.mVYY.Write("DV_DV_rawVariance");
  linReg.mVYYp.Write("DV_DV_rawVariance_prime");
  TVectorD mVY2(TMatrixDDiag(linReg.mVYY));
  mVY2.Write("DV_rawVariance");
  TVectorD mVP2(TMatrixDDiag(linReg.mVPP));
  mVP2.Write("IV_rawVariance");
  TVectorD mVY2prime(TMatrixDDiag(linReg.mVYYp));
  mVY2prime.Write("DV_rawVariance_prime");
 
  // normalized covariances
  linReg.mSPP.Write("IV_IV_normVariance");
  linReg.mSPY.Write("IV_DV_normVariance");
  linReg.mSYY.Write("DV_DV_normVariance");
  linReg.mSYYp.Write("DV_DV_normVariance_prime");
  TVectorD sigY2(TMatrixDDiag(linReg.mSYY));
  sigY2.Write("DV_normVariance");
  TVectorD sigX2(TMatrixDDiag(linReg.mSPP));
  sigX2.Write("IV_normVariance");
  TVectorD sigY2prime(TMatrixDDiag(linReg.mSYYp));
  sigY2prime.Write("DV_normVariance_prime");
 
  linReg.Axy.Write("A_xy");
  linReg.Ayx.Write("A_yx");
}