ATLAS Offline Software
Loading...
Searching...
No Matches
columnar::ColumnarPhysLiteTest Struct Reference

#include <ColumnarPhysliteTest.h>

Inheritance diagram for columnar::ColumnarPhysLiteTest:
Collaboration diagram for columnar::ColumnarPhysLiteTest:

Public Member Functions

 ColumnarPhysLiteTest ()
 ~ColumnarPhysLiteTest ()
std::string makeUniqueName ()
 make a unique tool name to be used in unit tests
void setupKnownColumns (std::span< const TestUtils::TestDefinition > testDefinitions)
void setupColumns (const ColumnVectorHeader &columnHeader)
void doCall (const TestUtils::TestDefinition &testDefinition)
void doCallMulti (const std::vector< TestUtils::TestDefinition > &testDefinitions)

Static Public Member Functions

static bool checkMode ()
 check whether we have the right mode

Public Attributes

std::unique_ptr< TFile > file
TTree * tree = nullptr
std::vector< std::shared_ptr< TestUtils::IColumnData > > knownColumns
std::vector< std::shared_ptr< TestUtils::IColumnData > > usedColumns
std::unordered_map< std::string, const std::vector< ColumnarOffsetType > * > offsetColumns

Detailed Description

Definition at line 38 of file ColumnarPhysliteTest.h.

Constructor & Destructor Documentation

◆ ColumnarPhysLiteTest()

columnar::ColumnarPhysLiteTest::ColumnarPhysLiteTest ( )

Definition at line 1794 of file ColumnarPhysliteTest.cxx.

1796 {
1797 static std::once_flag flag;
1798 std::call_once (flag, [] ()
1799 {
1800#ifdef XAOD_STANDALONE
1801 xAOD::Init().ignore();
1802#else
1803 POOL::Init();
1804#endif
1805 });
1806
1807 auto *fileName = getenv ("ASG_TEST_FILE_LITE_MC");
1808 if (fileName == nullptr)
1809 throw std::runtime_error ("missing ASG_TEST_FILE_LITE_MC");
1810 file.reset (TFile::Open (fileName, "READ"));
1811 if (!file)
1812 throw std::runtime_error ("failed to open file");
1813 tree = dynamic_cast<TTree*> (file->Get ("CollectionTree"));
1814 if (!tree)
1815 throw std::runtime_error ("failed to open tree");
1816 }
IAppMgrUI * Init(const char *options="POOLRootAccess/basic.opts")
Bootstraps (creates and configures) the Gaudi Application with the provided options file.
std::string getenv(const std::string &variableName)
get an environment variable
bool flag
Definition master.py:29
StatusCode Init(const char *appname)
Function initialising ROOT/PyROOT for using the ATLAS EDM.
Definition Init.cxx:31

◆ ~ColumnarPhysLiteTest()

columnar::ColumnarPhysLiteTest::~ColumnarPhysLiteTest ( )
default

Member Function Documentation

◆ checkMode()

bool columnar::ColumnarPhysLiteTest::checkMode ( )
static

check whether we have the right mode

Definition at line 1826 of file ColumnarPhysliteTest.cxx.

1828 {
1829 return true;
1830 }

◆ doCall()

void columnar::ColumnarPhysLiteTest::doCall ( const TestUtils::TestDefinition & testDefinition)

Definition at line 2073 of file ColumnarPhysliteTest.cxx.

2074 {
2075 doCallMulti ({testDefinition});
2076 }
void doCallMulti(const std::vector< TestUtils::TestDefinition > &testDefinitions)

◆ doCallMulti()

void columnar::ColumnarPhysLiteTest::doCallMulti ( const std::vector< TestUtils::TestDefinition > & testDefinitions)

Definition at line 2078 of file ColumnarPhysliteTest.cxx.

2079 {
2080 using namespace asg::msgUserCode;
2081 auto userConfiguration = TestUtils::UserConfiguration::fromEnvironment();
2082
2083 // apply systematics for all test definitions
2084 for (const auto& td : testDefinitions)
2085 {
2086 if (!td.sysName.empty())
2087 {
2088 auto *sysTool = dynamic_cast<CP::ISystematicsTool*>(td.tool);
2089 if (!sysTool)
2090 throw std::runtime_error ("tool does not support systematics");
2091 std::cout << "applying systematic variation: " << td.sysName << std::endl;
2092 if (sysTool->applySystematicVariation (CP::SystematicSet (td.sysName)).isFailure())
2093 throw std::runtime_error ("failed to apply systematic variation: " + td.sysName);
2094 }
2095 }
2096
2097 if constexpr (columnarAccessMode == 2)
2098 {
2099 // Create shared column header for all tools
2100 ColumnVectorHeader columnHeader;
2101
2102 // Build vector of ToolData from all testDefinitions
2103 std::vector<TestUtils::ToolData> toolDataVec;
2104 for (const auto& td : testDefinitions)
2105 toolDataVec.emplace_back (userConfiguration, td, columnHeader);
2106
2107 setupKnownColumns (testDefinitions);
2108 // Set up columns using the shared header (all tools have already registered
2109 // their columns via ToolColumnVectorMap, so we get all columns from the header)
2110 setupColumns (columnHeader);
2111
2112 // Connect column indices from header to each column for direct setting
2113 for (auto& column : usedColumns)
2114 column->connectColumnIndices (columnHeader);
2115
2116 Benchmark benchmarkEmpty ("empty");
2117 Benchmark benchmarkCheck ("", userConfiguration.batchSize);
2118
2119 const auto numberOfEvents = tree->GetEntries();
2120 Long64_t entry = 0;
2121 const auto startTime = std::chrono::high_resolution_clock::now();
2122 bool endLoop = false;
2123 for (; !endLoop; ++entry)
2124 {
2125 // just sample how much overhead there is for starting and
2126 // stopping the timer
2127 benchmarkEmpty.startTimer ();
2128 benchmarkEmpty.stopTimer ();
2129
2130 ColumnVectorData columnData (&columnHeader);
2131 for (auto& column : usedColumns)
2132 column->getEntry (entry % numberOfEvents);
2133 if ((entry + 1) % userConfiguration.batchSize == 0)
2134 {
2135 if (entry < numberOfEvents)
2136 {
2137 for (auto& column : usedColumns)
2138 column->collectColumnData ();
2139 }
2140 for (auto& column : usedColumns)
2141 column->setData (columnData);
2142 // Check data once (shared column data)
2143 benchmarkCheck.startTimer ();
2144 columnData.checkData ();
2145 benchmarkCheck.stopTimer ();
2146 // Call each tool
2147 for (auto& toolData : toolDataVec)
2148 toolData.call (columnData);
2149 for (auto& column : usedColumns)
2150 column->clearColumns ();
2151 if ((std::chrono::high_resolution_clock::now() - startTime) > userConfiguration.targetTime)
2152 endLoop = true;
2153 } else if (entry + 1 == numberOfEvents)
2154 {
2155 for (auto& column : usedColumns)
2156 column->collectColumnData ();
2157 }
2158 }
2159 std::cout << "Entries in file: " << numberOfEvents << std::endl;
2160 std::cout << "Total entries read: " << entry << std::endl;
2161 const float emptyTime = benchmarkEmpty.getEntryTime(0).value();
2162 std::cout << "Empty benchmark time: " << emptyTime << "ns (tick=" << Benchmark::getTickDuration() << "ns)" << std::endl;
2163 benchmarkEmpty.setSilence();
2164 const auto checkTime = benchmarkCheck.getEntryTime(emptyTime);
2165 if (checkTime)
2166 std::cout << "Check data time: " << checkTime.value() << "ns" << std::endl;
2167 benchmarkCheck.setSilence();
2168 {
2169 std::vector<TestUtils::BranchPerfData> branchPerfData;
2170 TestUtils::BranchPerfData summary;
2171 summary.name = "total";
2172 summary.timeRead = 0;
2173 summary.timeUnpack = 0;
2174 summary.timeShallowCopy = 0;
2175 summary.entries = std::nullopt;
2176 summary.nullEntries = std::nullopt;
2177 for (auto& column : usedColumns)
2178 {
2179 branchPerfData.push_back (column->getPerfData (emptyTime));
2180 summary.timeRead.value() += branchPerfData.back().timeRead.value_or(0);
2181 summary.timeUnpack.value() += branchPerfData.back().timeUnpack.value_or(0);
2182 summary.timeShallowCopy.value() += branchPerfData.back().timeShallowCopy.value_or(0);
2183 }
2184 std::sort (branchPerfData.begin(), branchPerfData.end(), [] (const auto& a, const auto& b) {return a.name < b.name;});
2185 branchPerfData.insert (branchPerfData.end(), summary);
2186 const std::size_t nameWidth = std::max_element (branchPerfData.begin(), branchPerfData.end(), [] (const auto& a, const auto& b) {return a.name.size() < b.name.size();})->name.size();
2187 std::string header = std::format ("{:{}} | read(ns) | unpack(ns) | size(B) | rate(MB/s) | compression | baskets | entries | null", "branch name", nameWidth);
2188 std::cout << "\n" << header << std::endl;
2189 std::cout << std::string (header.size(), '-') << std::endl;
2190 for (auto& data : branchPerfData)
2191 {
2192 if (data.name == "total")
2193 std::cout << std::string (header.size(), '-') << std::endl;
2194 std::cout << std::format ("{:{}} |", data.name, nameWidth);
2195 if (data.timeRead)
2196 std::cout << std::format ("{:>9.0f} |", data.timeRead.value());
2197 else
2198 std::cout << " |";
2199 if (data.timeUnpack)
2200 std::cout << std::format ("{:>11.1f} |", data.timeUnpack.value());
2201 else
2202 std::cout << " |";
2203 if (data.entrySize)
2204 std::cout << std::format ("{:>8.1f} |", data.entrySize.value());
2205 else
2206 std::cout << " |";
2207 if (data.timeRead && data.entrySize)
2208 std::cout << std::format ("{:>11.1f} |", (data.entrySize.value() / (data.timeRead.value() * 1e-3 * 1.024 * 1.024)));
2209 else
2210 std::cout << " |";
2211 if (data.entrySize && data.uncompressedSize)
2212 std::cout << std::format ("{:>12.2f} |", float (data.uncompressedSize.value()) / data.entrySize.value());
2213 else
2214 std::cout << " |";
2215 if (data.numBaskets)
2216 std::cout << std::format ("{:>8} |", data.numBaskets.value());
2217 else
2218 std::cout << " |";
2219 if (data.entries)
2220 std::cout << std::format ("{:>8.2f} |", static_cast<float>(data.entries.value())/numberOfEvents);
2221 else
2222 std::cout << " |";
2223 if (data.nullEntries && data.entries)
2224 std::cout << std::format ("{:>4.0f}%", static_cast<float>(data.nullEntries.value()) / data.entries.value() * 100.0f);
2225 std::cout << std::endl;
2226 }
2227 }
2228 {
2229 std::vector<TestUtils::ToolPerfData> toolPerfData;
2230 for (auto& toolData : toolDataVec)
2231 {
2232 toolPerfData.emplace_back ();
2233 toolPerfData.back().name = toolData.name;
2234 toolPerfData.back().timeCall = toolData.benchmarkCall.getEntryTime (emptyTime);
2235 if (userConfiguration.runToolTwice)
2236 toolPerfData.back().timeCall2 = toolData.benchmarkCall2.getEntryTime (emptyTime);
2237 }
2238 const std::size_t nameWidth = std::max_element (toolPerfData.begin(), toolPerfData.end(), [] (const auto& a, const auto& b) {return a.name.size() < b.name.size();})->name.size();
2239 std::string header = std::format ("{:{}} | call(ns) | call2(ns)", "tool name", nameWidth);
2240 std::cout << "\n" << header << std::endl;
2241 std::cout << std::string (header.size(), '-') << std::endl;
2242 for (auto& data : toolPerfData)
2243 {
2244 std::cout << std::format ("{:{}} |", data.name, nameWidth);
2245 if (data.timeCall)
2246 std::cout << std::format ("{:>9.0f} |", data.timeCall.value());
2247 else
2248 std::cout << " |";
2249 if (data.timeCall2)
2250 std::cout << std::format ("{:>10.0f}", data.timeCall2.value());
2251 else
2252 std::cout << " ";
2253 std::cout << std::endl;
2254 }
2255 // Add totals line for multiple tools
2256 if (toolPerfData.size() > 1)
2257 {
2258 std::optional<float> totalCall, totalCall2;
2259 for (const auto& data : toolPerfData)
2260 {
2261 if (data.timeCall)
2262 totalCall = totalCall.value_or (0) + data.timeCall.value();
2263 if (data.timeCall2)
2264 totalCall2 = totalCall2.value_or (0) + data.timeCall2.value();
2265 }
2266 std::cout << std::string (header.size(), '-') << std::endl;
2267 std::cout << std::format ("{:{}} |", "total", nameWidth);
2268 if (totalCall)
2269 std::cout << std::format ("{:>9.0f} |", totalCall.value());
2270 else
2271 std::cout << " |";
2272 if (totalCall2)
2273 std::cout << std::format ("{:>10.0f}", totalCall2.value());
2274 else
2275 std::cout << " ";
2276 std::cout << std::endl;
2277 }
2278 }
2279 } else if constexpr (columnarAccessMode == 0)
2280 {
2281 TestUtils::runXaodTest (userConfiguration, testDefinitions, file.get());
2282 } else if constexpr (columnarAccessMode == 100)
2283 {
2284 const auto& testDefinition = testDefinitions[0];
2285 TestUtils::runXaodArrayTest (userConfiguration, testDefinition, file.get());
2286 }
2287 }
void checkTime()
int numberOfEvents()
char data[hepevt_bytes_allocation_ATLAS]
Definition HepEvt.cxx:11
static Double_t a
static float getTickDuration()
Definition Benchmark.h:86
void runXaodArrayTest(const UserConfiguration &userConfiguration, const TestDefinition &testDefinition, TFile *file)
void runXaodTest(const UserConfiguration &userConfiguration, std::span< const TestDefinition > testDefinitions, TFile *file)
constexpr unsigned columnarAccessMode
Definition ColumnarDef.h:15
void sort(typename DataModel_detail::iterator< DVL > beg, typename DataModel_detail::iterator< DVL > end)
Specialization of sort for DataVector/List.
std::vector< std::shared_ptr< TestUtils::IColumnData > > usedColumns
void setupKnownColumns(std::span< const TestUtils::TestDefinition > testDefinitions)
void setupColumns(const ColumnVectorHeader &columnHeader)
static UserConfiguration fromEnvironment()
create a UserConfiguration, loading from the file pointed to by the COLUMNAR_TEST_CONFIG environment ...

◆ makeUniqueName()

std::string columnar::ColumnarPhysLiteTest::makeUniqueName ( )

make a unique tool name to be used in unit tests

Definition at line 1820 of file ColumnarPhysliteTest.cxx.

1821 {
1822 static std::atomic<unsigned> index = 0;
1823 return "UniquePhysliteTestTool" + std::to_string(++index);
1824 }
str index
Definition DeMoScan.py:362

◆ setupColumns()

void columnar::ColumnarPhysLiteTest::setupColumns ( const ColumnVectorHeader & columnHeader)

Definition at line 1999 of file ColumnarPhysliteTest.cxx.

2000 {
2001 using namespace asg::msgUserCode;
2002
2003 // Get all column info directly from the header (all tools have already
2004 // registered their columns via ToolColumnVectorMap)
2005 auto requestedColumns = columnHeader.getAllColumnInfo();
2006
2007 // Print requested columns
2008 for (auto& [name, info] : requestedColumns)
2009 std::cout << "requested columns: " << name << std::endl;
2010
2011 for (auto& column : knownColumns)
2012 {
2013 if (column->connect (tree, offsetColumns, requestedColumns))
2014 usedColumns.push_back (column);
2015 }
2016
2017 std::set<std::string> unclaimedColumns;
2018 for (auto& column : requestedColumns)
2019 {
2020 if (!column.second.isOptional)
2021 unclaimedColumns.insert (column.first);
2022 else
2023 std::cout << "optional column not claimed: " << column.first << std::endl;
2024 }
2025 std::erase_if (unclaimedColumns, [&] (auto& columnName)
2026 {
2027 const auto& info = requestedColumns.at (columnName);
2028 if (info.accessMode != ColumnAccessMode::output || !info.fixedDimensions.empty())
2029 return false;
2030 auto offsetIter = std::find_if (usedColumns.begin(), usedColumns.end(), [&] (const std::shared_ptr<TestUtils::IColumnData>& column)
2031 {
2032 for (auto& output : column->outputColumns)
2033 {
2034 if (output.name == info.offsetName)
2035 return true;
2036 }
2037 return false;
2038 });
2039 if (offsetIter == usedColumns.end())
2040 return false;
2041 std::shared_ptr<TestUtils::IColumnData> myColumn;
2042 if (*info.type == typeid(float))
2043 myColumn = std::make_shared<TestUtils::ColumnDataOutVector<float>> (info.name, 0);
2044 else if (*info.type == typeid(char))
2045 myColumn = std::make_shared<TestUtils::ColumnDataOutVector<char>> (info.name, 0);
2046 else if (*info.type == typeid(std::uint16_t))
2047 myColumn = std::make_shared<TestUtils::ColumnDataOutVector<std::uint16_t>> (info.name, 0);
2048 else if (*info.type == typeid(std::uint64_t))
2049 myColumn = std::make_shared<TestUtils::ColumnDataOutVector<std::uint64_t>> (info.name, 0);
2050 else
2051 {
2052 ANA_MSG_WARNING ("unhandled column type: " << info.name << " " << info.type->name());
2053 return false;
2054 }
2055 knownColumns.push_back (myColumn);
2056 if (!myColumn->connect (tree, offsetColumns, requestedColumns))
2057 {
2058 ANA_MSG_WARNING ("failed to connect dynamic output column: " << info.name);
2059 return false;
2060 }
2061 usedColumns.push_back (myColumn);
2062 return true;
2063 });
2064 if (!unclaimedColumns.empty())
2065 {
2066 std::string message = "columns not claimed:";
2067 for (auto& column : unclaimedColumns)
2068 message += " " + column;
2069 throw std::runtime_error (message);
2070 }
2071 }
#define ANA_MSG_WARNING(xmsg)
Macro printing warning messages.
@ output
an output column
Definition ColumnInfo.h:24
std::size_t erase_if(T_container &container, T_Func pred)
std::vector< std::shared_ptr< TestUtils::IColumnData > > knownColumns
std::unordered_map< std::string, const std::vector< ColumnarOffsetType > * > offsetColumns

◆ setupKnownColumns()

void columnar::ColumnarPhysLiteTest::setupKnownColumns ( std::span< const TestUtils::TestDefinition > testDefinitions)

Definition at line 1832 of file ColumnarPhysliteTest.cxx.

1833 {
1834 using namespace TestUtils;
1835
1836 knownColumns.push_back (std::make_shared<ColumnDataEventCount> ());
1837
1838 tree->SetMakeClass (1);
1839 {
1840 std::unordered_map<std::string,TBranch*> branches;
1841 {
1842 TIter branchIter (tree->GetListOfBranches());
1843 TObject *obj = nullptr;
1844 while ((obj = branchIter()))
1845 {
1846 TBranch *branch = nullptr;
1847 if ((branch = dynamic_cast<TBranch*>(obj)))
1848 {
1849 branches.emplace (branch->GetName(), branch);
1850 TIter subBranchIter (branch->GetListOfBranches());
1851 while ((obj = subBranchIter()))
1852 {
1853 if (auto subBranch = dynamic_cast<TBranch*>(obj))
1854 branches.emplace (subBranch->GetName(), subBranch);
1855 }
1856 }
1857 }
1858 }
1859
1860 for (const auto& [name, branch] : branches)
1861 {
1862 if (name.find ("AuxDyn.") != std::string::npos ||
1863 name.find ("Aux.") != std::string::npos)
1864 {
1865 TClass *branchClass = nullptr;
1866 EDataType branchType {};
1867 branch->GetExpectedType (branchClass, branchType);
1868 if (branchClass == nullptr)
1869 {
1870 switch (branchType)
1871 {
1872 case kInt_t:
1873 knownColumns.push_back (std::make_shared<ColumnDataScalar<std::int32_t>> (branch->GetName()));
1874 break;
1875 case kUInt_t:
1876 knownColumns.push_back (std::make_shared<ColumnDataScalar<std::uint32_t>> (branch->GetName()));
1877 break;
1878 case kULong_t:
1879 knownColumns.push_back (std::make_shared<ColumnDataScalar<std::uint64_t>> (branch->GetName()));
1880 break;
1881 case kULong64_t:
1882 knownColumns.push_back (std::make_shared<ColumnDataScalar<std::uint64_t>> (branch->GetName()));
1883 break;
1884 case kFloat_t:
1885 knownColumns.push_back (std::make_shared<ColumnDataScalar<float>> (branch->GetName()));
1886 break;
1887 default:
1888 // no-op
1889 break;
1890 }
1891 } else
1892 {
1893 if (*branchClass->GetTypeInfo() == typeid(std::vector<float>))
1894 {
1895 knownColumns.push_back (std::make_shared<ColumnDataVector<float>> (branch->GetName()));
1896 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<char>))
1897 {
1898 knownColumns.push_back (std::make_shared<ColumnDataVector<char>> (branch->GetName()));
1899 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::int8_t>))
1900 {
1901 knownColumns.push_back (std::make_shared<ColumnDataVector<std::int8_t>> (branch->GetName()));
1902 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::uint8_t>))
1903 {
1904 knownColumns.push_back (std::make_shared<ColumnDataVector<std::uint8_t>> (branch->GetName()));
1905 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::int16_t>))
1906 {
1907 knownColumns.push_back (std::make_shared<ColumnDataVector<std::int16_t>> (branch->GetName()));
1908 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::uint16_t>))
1909 {
1910 knownColumns.push_back (std::make_shared<ColumnDataVector<std::uint16_t>> (branch->GetName()));
1911 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::int32_t>))
1912 {
1913 knownColumns.push_back (std::make_shared<ColumnDataVector<std::int32_t>> (branch->GetName()));
1914 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::uint32_t>))
1915 {
1916 knownColumns.push_back (std::make_shared<ColumnDataVector<std::uint32_t>> (branch->GetName()));
1917 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::int64_t>))
1918 {
1919 knownColumns.push_back (std::make_shared<ColumnDataVector<std::int64_t>> (branch->GetName()));
1920 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::uint64_t>))
1921 {
1922 knownColumns.push_back (std::make_shared<ColumnDataVector<std::uint64_t>> (branch->GetName()));
1923 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::vector<float>>))
1924 {
1925 knownColumns.push_back (std::make_shared<ColumnDataVectorVector<float>> (branch->GetName()));
1926 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::vector<std::int32_t>>))
1927 {
1928 knownColumns.push_back (std::make_shared<ColumnDataVectorVector<std::int32_t>> (branch->GetName()));
1929 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::vector<std::uint64_t>>))
1930 {
1931 knownColumns.push_back (std::make_shared<ColumnDataVectorVector<std::uint64_t>> (branch->GetName()));
1932 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::vector<std::vector<std::size_t>>>))
1933 {
1934 knownColumns.push_back (std::make_shared<ColumnDataVectorVectorVector<std::size_t>> (branch->GetName()));
1935 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::vector<std::vector<unsigned char>>>))
1936 {
1937 knownColumns.push_back (std::make_shared<ColumnDataVectorVectorVector<unsigned char>> (branch->GetName()));
1938 } else if (*branchClass->GetTypeInfo() == typeid(std::vector<std::string>))
1939 {
1940 knownColumns.push_back (std::make_shared<ColumnDataMetNames> (branch->GetName()));
1941 }
1942 }
1943 }
1944 }
1945 }
1946
1947 // This is a fallback for the case that we don't have an explicit
1948 // `samplingPattern` branch in our input file (i.e. an older file),
1949 // to allow us to still test tools needing it. This is likely not
1950 // something that actual users can do (they need the new files), but
1951 // for testing it seems like a reasonable workaround.
1952 knownColumns.push_back (std::make_shared<ColumnDataSamplingPattern> ("egammaClusters"));
1953
1954 // For branches that are element links they need to be explicitly
1955 // declared to have the correct xAOD type, correct split setting,
1956 // and correct linked containers.
1957 knownColumns.push_back (std::make_shared<ColumnDataVectorVectorLink<xAOD::CaloClusterContainer>> ("AnalysisElectronsAuxDyn.caloClusterLinks"));
1958 knownColumns.push_back (std::make_shared<ColumnDataVectorVectorLink<xAOD::TrackParticleContainer>> ("AnalysisElectronsAuxDyn.trackParticleLinks"));
1959 knownColumns.push_back (std::make_shared<ColumnDataVectorVectorLink<xAOD::CaloClusterContainer>> ("AnalysisPhotonsAuxDyn.caloClusterLinks"));
1960 knownColumns.push_back (std::make_shared<ColumnDataVectorVectorLink<xAOD::VertexContainer>> ("AnalysisPhotonsAuxDyn.vertexLinks"));
1961 knownColumns.push_back (std::make_shared<ColumnDataVectorSplitLink<xAOD::TrackParticleContainer>> ("AnalysisMuonsAuxDyn.inDetTrackParticleLink"));
1962 knownColumns.push_back (std::make_shared<ColumnDataVectorSplitLink<xAOD::TrackParticleContainer>> ("AnalysisMuonsAuxDyn.combinedTrackParticleLink"));
1963 knownColumns.push_back (std::make_shared<ColumnDataVectorSplitLink<xAOD::TrackParticleContainer>> ("AnalysisMuonsAuxDyn.extrapolatedMuonSpectrometerTrackParticleLink"));
1964 knownColumns.push_back (std::make_shared<ColumnDataVectorVectorLink<xAOD::TrackParticleContainer>> ("GSFConversionVerticesAuxDyn.trackParticleLinks"));
1965 knownColumns.push_back (std::make_shared<ColumnDataVectorSplitLink<xAOD::TrackParticleContainer>> ("GSFTrackParticlesAuxDyn.originalTrackParticle"));
1966 knownColumns.push_back (std::make_shared<ColumnDataVectorVectorVariantLink<xAOD::IParticleContainer>>("AnalysisJetsAuxDyn.GhostTrack"));
1967 knownColumns.push_back (std::make_shared<ColumnDataVectorLink<xAOD::JetContainer>>("METAssoc_AnalysisMETAux.jetLink"));
1968 knownColumns.push_back (std::make_shared<ColumnDataVectorVectorVariantLink<xAOD::IParticleContainer>>("METAssoc_AnalysisMETAux.objectLinks"));
1969
1970 // For METMaker we need to preplace all of the MET terms that we
1971 // expect to be used, that's what this line does.
1972 std::vector<std::string> allMetTermNames;
1973 for (const auto& td : testDefinitions)
1974 {
1975 for (const auto& name : td.metTermNames)
1976 {
1977 if (std::find (allMetTermNames.begin(), allMetTermNames.end(), name) == allMetTermNames.end())
1978 allMetTermNames.push_back (name);
1979 }
1980 }
1981 if (!allMetTermNames.empty())
1982 knownColumns.push_back (std::make_shared<ColumnDataOutputMet> ("OutputMET", allMetTermNames));
1983
1984 // For METMaker we need various extra columns to run. This may need
1985 // some work to avoid, but would likey be worth it.
1986 knownColumns.push_back (std::make_shared<ColumnDataOutVector<std::uint16_t>> ("AnalysisMuons.objectType", xAOD::Type::Muon));
1987 knownColumns.push_back (std::make_shared<ColumnDataOutVector<float>> ("AnalysisMuons.m", ParticleConstants::muonMassInMeV));
1988 knownColumns.push_back (std::make_shared<ColumnDataOutVector<std::uint16_t>> ("AnalysisJets.objectType", xAOD::Type::Jet));
1989
1990 // These are columns that represent variables that are normally held
1991 // by METAssociationHelper, or alternatively are decorated on the
1992 // MET terms (even though they are per object).
1993 knownColumns.push_back (std::make_shared<ColumnDataOutVector<float>> ("AnalysisMuons.MetObjectWeight", 0));
1994 knownColumns.push_back (std::make_shared<ColumnDataOutVector<float>> ("AnalysisJets.MetObjectWeight", 0));
1995 knownColumns.push_back (std::make_shared<ColumnDataOutVector<float>> ("AnalysisJets.MetObjectWeightSoft", 0));
1996 knownColumns.push_back (std::make_shared<ColumnDataOutVector<MissingETBase::Types::bitmask_t>> ("METAssoc_AnalysisMET.useObjectFlags", 0));
1997 }
constexpr double muonMassInMeV
the mass of the muon (in MeV)
@ Jet
The object is a jet.
Definition ObjectType.h:40
@ Muon
The object is a muon.
Definition ObjectType.h:48

Member Data Documentation

◆ file

std::unique_ptr<TFile> columnar::ColumnarPhysLiteTest::file

Definition at line 40 of file ColumnarPhysliteTest.h.

◆ knownColumns

std::vector<std::shared_ptr<TestUtils::IColumnData> > columnar::ColumnarPhysLiteTest::knownColumns

Definition at line 43 of file ColumnarPhysliteTest.h.

◆ offsetColumns

std::unordered_map<std::string,const std::vector<ColumnarOffsetType>*> columnar::ColumnarPhysLiteTest::offsetColumns

Definition at line 45 of file ColumnarPhysliteTest.h.

◆ tree

TTree* columnar::ColumnarPhysLiteTest::tree = nullptr

Definition at line 41 of file ColumnarPhysliteTest.h.

◆ usedColumns

std::vector<std::shared_ptr<TestUtils::IColumnData> > columnar::ColumnarPhysLiteTest::usedColumns

Definition at line 44 of file ColumnarPhysliteTest.h.


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