File indexing completed on 2024-04-06 12:23:38
0001 #include <iostream>
0002 #include <sstream>
0003 #include <fstream>
0004 #include <algorithm>
0005 #include <iterator>
0006 #include <cstddef>
0007 #include <cstring>
0008 #include <vector>
0009 #include <set>
0010
0011
0012 #include <RVersion.h>
0013
0014 #include <TBufferFile.h>
0015 #include <TClass.h>
0016
0017 #include "FWCore/Utilities/interface/Exception.h"
0018
0019 #include "PhysicsTools/MVAComputer/interface/zstream.h"
0020
0021 #include "PhysicsTools/MVAComputer/interface/MVAComputer.h"
0022 #include "PhysicsTools/MVAComputer/interface/VarProcessor.h"
0023 #include "PhysicsTools/MVAComputer/interface/Calibration.h"
0024 #include "PhysicsTools/MVAComputer/interface/Variable.h"
0025 #include "PhysicsTools/MVAComputer/interface/AtomicId.h"
0026
0027
0028
0029 #ifdef DEBUG_EVAL
0030 #include "FWCore/Utilities/interface/TypeDemangler.h"
0031 #endif
0032
0033 #define STANDALONE_HEADER "MVAComputer calibration\n"
0034
0035 namespace PhysicsTools {
0036
0037 MVAComputer::MVAComputer(const Calibration::MVAComputer *calib) : nVars(0), output(0) { setup(calib); }
0038
0039 MVAComputer::MVAComputer(Calibration::MVAComputer *calib, bool owned) : nVars(0), output(0) {
0040 if (owned)
0041 this->owned.reset(calib);
0042 setup(calib);
0043 }
0044
0045 MVAComputer::MVAComputer(const char *filename) : nVars(0), output(0), owned(readCalibration(filename)) {
0046 setup(owned.get());
0047 }
0048
0049 MVAComputer::MVAComputer(std::istream &is) : nVars(0), output(0), owned(readCalibration(is)) { setup(owned.get()); }
0050
0051 void MVAComputer::setup(const Calibration::MVAComputer *calib) {
0052 nVars = calib->inputSet.size();
0053 output = calib->output;
0054
0055 std::vector<Variable::Flags> flags(nVars, Variable::FLAG_ALL);
0056 const TrainMVAComputerCalibration *trainCalib = dynamic_cast<const TrainMVAComputerCalibration *>(calib);
0057 if (trainCalib)
0058 trainCalib->initFlags(flags);
0059
0060 VarProcessor::ConfigCtx config(flags);
0061 std::vector<const Calibration::VarProcessor *> processors = calib->getProcessors();
0062
0063 for (const auto *calibProc : processors) {
0064 std::string name = calibProc->getInstanceName();
0065 VarProcessor *processor = VarProcessor::create(name.c_str(), calibProc, this);
0066 if (!processor)
0067 throw cms::Exception("UnknownProcessor") << name << " could not be instantiated." << std::endl;
0068
0069 VarProcessor::ConfigCtx::iterator::difference_type pos = config.end() - config.begin();
0070 processor->configure(config);
0071 unsigned int nOutput = (config.end() - config.begin()) - pos;
0072 if (!nOutput)
0073 throw cms::Exception("InvalidProcessor") << name
0074 << " rejected input variable "
0075 "configuration"
0076 << std::endl;
0077
0078 varProcessors.push_back(Processor(processor, nOutput));
0079 }
0080
0081 for (VarProcessor::ConfigCtx::iterator iter = config.begin() + nVars; iter != config.end(); iter++) {
0082 VarProcessor::Config *origin = &config[iter->origin];
0083 if (iter->origin >= nVars)
0084 iter->origin = origin->origin;
0085
0086 if (iter->mask & Variable::FLAG_MULTIPLE) {
0087 iter->mask = (Variable::Flags)(iter->mask & origin->mask);
0088 config[iter->origin].origin++;
0089 }
0090 }
0091
0092 nVars = config.size();
0093
0094 if (output >= nVars)
0095
0096 throw cms::Exception("InvalidOutput") << "Output variable at index " << output << " invalid." << std::endl;
0097
0098 std::set<InputVar> variables;
0099 unsigned int i = 0;
0100 for (std::vector<Calibration::Variable>::const_iterator iter = calib->inputSet.begin();
0101 iter != calib->inputSet.end();
0102 ++iter, i++) {
0103 InputVar var;
0104 var.var = Variable(iter->name, config[i].mask);
0105 var.index = i;
0106 var.multiplicity = 0;
0107 variables.insert(var);
0108 }
0109
0110 inputVariables.resize(i);
0111
0112 std::copy(variables.begin(), variables.end(), inputVariables.begin());
0113
0114 for (unsigned int j = 0; j < i; j++)
0115 inputVariables[j].multiplicity = config[j].origin;
0116 }
0117
0118 MVAComputer::~MVAComputer() {}
0119
0120 int MVAComputer::getVariableId(AtomicId name) const {
0121 std::vector<InputVar>::const_iterator pos = std::lower_bound(inputVariables.begin(), inputVariables.end(), name);
0122
0123 if (pos == inputVariables.end() || pos->var.getName() != name)
0124 return -1;
0125 else
0126 return pos->index;
0127 }
0128
0129 template <class T>
0130 void MVAComputer::evalInternal(T &ctx) const {
0131 double *output = ctx.values() + ctx.n();
0132 int *outConf = ctx.conf() + inputVariables.size();
0133
0134 #ifdef DEBUG_EVAL
0135 std::cout << "Input" << std::endl;
0136 double *v = ctx.values();
0137 for (int *o = ctx.conf(); o < outConf; o++) {
0138 std::cout << "\tVar " << (o - ctx.conf()) << std::endl;
0139 for (int i = o[0]; i < o[1]; i++)
0140 std::cout << "\t\t" << *v++ << std::endl;
0141 }
0142 #endif
0143 std::vector<Processor>::const_iterator iter = varProcessors.begin();
0144 while (iter != varProcessors.end()) {
0145 std::vector<Processor>::const_iterator loop = iter;
0146 int *loopOutConf = outConf;
0147 int *loopStart = nullptr;
0148 double *loopOutput = output;
0149 VarProcessor::LoopCtx loopCtx;
0150
0151 VarProcessor::LoopStatus status = VarProcessor::kNext;
0152 unsigned int offset = 0;
0153 while (status != VarProcessor::kStop) {
0154 std::vector<Processor>::const_iterator next = iter + 1;
0155 unsigned int nextOutput = (next != varProcessors.end()) ? next->nOutput : 0;
0156
0157 #ifdef DEBUG_EVAL
0158 std::string demangledName;
0159 edm::typeDemangle(typeid(*iter->processor).name(), demangledName);
0160 std::cout << demangledName << std::endl;
0161 #endif
0162 if (status != VarProcessor::kSkip)
0163 ctx.eval(
0164 &*iter->processor, outConf, output, loopStart ? loopStart : loopOutConf, loopCtx, offset, iter->nOutput);
0165
0166 #ifdef DEBUG_EVAL
0167 for (unsigned int i = 0; i < iter->nOutput; i++, outConf++) {
0168 std::cout << "\tVar " << (outConf - ctx.conf()) << std::endl;
0169 for (int j = outConf[0]; j < outConf[1]; j++)
0170 std::cout << "\t\t" << *output++ << std::endl;
0171 }
0172 #else
0173 int orig = *outConf;
0174 outConf += iter->nOutput;
0175 output += *outConf - orig;
0176 #endif
0177
0178 status = loop->processor->loop(output, outConf, nextOutput, loopCtx, offset);
0179
0180 if (status == VarProcessor::kReset) {
0181 outConf = loopOutConf;
0182 output = loopOutput;
0183 loopStart = nullptr;
0184 offset = 0;
0185 iter = loop;
0186 } else {
0187 if (loop == iter)
0188 loopStart = outConf;
0189 iter = next;
0190 }
0191 }
0192 }
0193 }
0194
0195 Calibration::MVAComputer *MVAComputer::readCalibration(const char *filename) {
0196 std::ifstream file(filename);
0197 return readCalibration(file);
0198 }
0199
0200 Calibration::MVAComputer *MVAComputer::readCalibration(std::istream &is) {
0201 if (!is.good())
0202 throw cms::Exception("InvalidFileState") << "Stream passed to MVAComputer::readCalibration "
0203 "has an invalid state."
0204 << std::endl;
0205
0206 char header[sizeof STANDALONE_HEADER - 1] = {
0207 0,
0208 };
0209 if (is.readsome(header, sizeof header) != sizeof header ||
0210 std::memcmp(header, STANDALONE_HEADER, sizeof header) != 0)
0211 throw cms::Exception("InvalidFileFormat") << "Stream passed to MVAComputer::readCalibration "
0212 "is not a valid calibration file."
0213 << std::endl;
0214
0215 TClass *rootClass = TClass::GetClass("PhysicsTools::Calibration::MVAComputer");
0216 if (!rootClass)
0217 throw cms::Exception("DictionaryMissing") << "CondFormats dictionary for "
0218 "PhysicsTools::Calibration::MVAComputer missing"
0219 << std::endl;
0220
0221 ext::izstream izs(&is);
0222 std::ostringstream ss;
0223 ss << izs.rdbuf();
0224 std::string buf = ss.str();
0225
0226 TBufferFile buffer(TBuffer::kRead, buf.size(), const_cast<void *>(static_cast<const void *>(buf.c_str())), kFALSE);
0227 buffer.InitMap();
0228
0229 std::unique_ptr<Calibration::MVAComputer> calib(new Calibration::MVAComputer());
0230 buffer.StreamObject(static_cast<void *>(calib.get()), rootClass);
0231
0232 return calib.release();
0233 }
0234
0235 void MVAComputer::writeCalibration(const char *filename, const Calibration::MVAComputer *calib) {
0236 std::ofstream file(filename);
0237 writeCalibration(file, calib);
0238 }
0239
0240 void MVAComputer::writeCalibration(std::ostream &os, const Calibration::MVAComputer *calib) {
0241 if (!os.good())
0242 throw cms::Exception("InvalidFileState") << "Stream passed to MVAComputer::writeCalibration "
0243 "has an invalid state."
0244 << std::endl;
0245
0246 os << STANDALONE_HEADER;
0247
0248 TClass *rootClass = TClass::GetClass("PhysicsTools::Calibration::MVAComputer");
0249 if (!rootClass)
0250 throw cms::Exception("DictionaryMissing") << "CondFormats dictionary for "
0251 "PhysicsTools::Calibration::MVAComputer missing"
0252 << std::endl;
0253
0254 TBufferFile buffer(TBuffer::kWrite);
0255 buffer.StreamObject(const_cast<void *>(static_cast<const void *>(calib)), rootClass);
0256
0257 ext::ozstream ozs(&os);
0258 ozs.write(buffer.Buffer(), buffer.Length());
0259 ozs.flush();
0260 }
0261
0262
0263
0264 void MVAComputer::DerivContext::eval(const VarProcessor *proc,
0265 int *outConf,
0266 double *output,
0267 int *loop,
0268 VarProcessor::LoopCtx &ctx,
0269 unsigned int offset,
0270 unsigned int out) const {
0271 proc->deriv(values(), conf(), output, outConf, loop, ctx, offset, n(), out, deriv_);
0272 }
0273
0274 double MVAComputer::DerivContext::output(unsigned int output, std::vector<double> &derivs) const {
0275 derivs.clear();
0276 derivs.reserve(n_);
0277
0278 int pos = conf_[output];
0279 if (pos >= (int)n_)
0280 std::copy(&deriv_.front() + n_ * (pos - n_), &deriv_.front() + n_ * (pos - n_ + 1), std::back_inserter(derivs));
0281 else {
0282 derivs.resize(n_);
0283 derivs[pos] = 1.;
0284 }
0285
0286 return values_[pos];
0287 }
0288
0289 template void MVAComputer::evalInternal(EvalContext &ctx) const;
0290 template void MVAComputer::evalInternal(DerivContext &ctx) const;
0291
0292 }
