Line Code
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120
#include "CalibTracker/SiStripChannelGain/plugins/SiStripGainFromAsciiFile.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"

#include "CalibTracker/SiStripCommon/interface/SiStripDetInfoFileReader.h"

#include "DataFormats/DetId/interface/DetId.h"

#include <iostream>
#include <fstream>
#include <sstream>

SiStripGainFromAsciiFile::SiStripGainFromAsciiFile(const edm::ParameterSet& iConfig)
    : ConditionDBWriter<SiStripApvGain>(iConfig) {
  Asciifilename_ = iConfig.getParameter<std::string>("InputFileName");
  referenceValue_ = iConfig.getParameter<double>("referenceValue");
  fp_ = iConfig.getUntrackedParameter<edm::FileInPath>("file", edm::FileInPath(SiStripDetInfoFileReader::kDefaultFile));
}

SiStripGainFromAsciiFile::~SiStripGainFromAsciiFile() {
  edm::LogInfo("SiStripGainFromAsciiFile::~SiStripGainFromAsciiFile");
}

std::unique_ptr<SiStripApvGain> SiStripGainFromAsciiFile::getNewObject() {
  auto obj = std::make_unique<SiStripApvGain>();

  std::stringstream ss;
  ss.str("");
  ss << "[SiStripGainFromAsciiFile::getNewObject]\n Reading Ascii File\n";
  FILE* infile = fopen(Asciifilename_.c_str(), "r");
  char line[4096];
  if (infile) {
    while (fgets(line, 4096, infile) != nullptr) {
      uint32_t detid;
      ModuleGain MG;
      MG.apv[0] = 0.0;
      MG.apv[1] = 0.0;
      MG.apv[2] = 0.0;
      MG.apv[3] = 0.0;
      MG.apv[4] = 0.0;
      MG.apv[5] = 0.0;
      char* saveptr;
      char* pch = strtok_r(line, " ", &saveptr);
      int Arg = 0;
      while (pch != nullptr) {
        if (Arg == 0) {
          sscanf(pch, "%d", &detid);
        } else if (Arg <= 6) {
          sscanf(pch, "%f", &(MG.apv[Arg - 1]));
        } else {
          //nothing to do here
        }
        pch = strtok_r(nullptr, " ", &saveptr);
        Arg++;
      }
      ss << detid << " " << MG.apv[0] << " " << MG.apv[1] << " " << MG.apv[2] << " " << MG.apv[3] << " " << MG.apv[4]
         << " " << MG.apv[5] << std::endl;
      GainsMap.insert(std::pair<unsigned int, ModuleGain>(detid, MG));
    }
    fclose(infile);
    edm::LogInfo("SiStripGainFromAsciiFile") << ss.str();
  } else {
    edm::LogError("SiStripGainFromAsciiFile")
        << " [SiStripGainFromAsciiFile::getNewObject] Error opening file " << Asciifilename_ << std::endl;
    assert(0);
  }

  const auto detInfo = SiStripDetInfoFileReader::read(fp_.fullPath());

  ss.str("");
  ss << "[SiStripGainFromAsciiFile::getNewObject]\n Filling SiStripApvGain object";
  short nApvPair;
  for (const auto it : detInfo.getAllDetIds()) {
    ModuleGain MG;
    if (DetId(it).det() != DetId::Tracker)
      continue;

    nApvPair = detInfo.getNumberOfApvsAndStripLength(it).first / 2;

    ss << "Looking at detid " << it << " nApvPairs  " << nApvPair << std::endl;
    auto iter = GainsMap.find(it);
    if (iter != GainsMap.end()) {
      MG = iter->second;
      ss << " " << MG.apv[0] << " " << MG.apv[1] << " " << MG.apv[2] << " " << MG.apv[3] << " " << MG.apv[4] << " "
         << MG.apv[5] << std::endl;
    } else {
      ss << "Hard reset for detid " << it << std::endl;
      MG.hard_reset(referenceValue_);
    }

    std::vector<float> DetGainsVector;

    if (nApvPair == 2) {
      DetGainsVector.push_back(MG.apv[0] / referenceValue_);
      DetGainsVector.push_back(MG.apv[1] / referenceValue_);
      DetGainsVector.push_back(MG.apv[2] / referenceValue_);
      DetGainsVector.push_back(MG.apv[3] / referenceValue_);
    } else if (nApvPair == 3) {
      DetGainsVector.push_back(MG.apv[0] / referenceValue_);
      DetGainsVector.push_back(MG.apv[1] / referenceValue_);
      DetGainsVector.push_back(MG.apv[2] / referenceValue_);
      DetGainsVector.push_back(MG.apv[3] / referenceValue_);
      DetGainsVector.push_back(MG.apv[4] / referenceValue_);
      DetGainsVector.push_back(MG.apv[5] / referenceValue_);
    } else {
      edm::LogError("SiStripGainFromAsciiFile") << " SiStripGainFromAsciiFile::getNewObject] ERROR for detid " << it
                                                << " not expected number of APV pairs " << nApvPair << std::endl;
    }

    SiStripApvGain::Range range(DetGainsVector.begin(), DetGainsVector.end());
    if (!obj->put(it, range)) {
      edm::LogError("SiStripGainFromAsciiFile")
          << " [SiStripGainFromAsciiFile::getNewObject] detid already exists" << std::endl;
      ss << " [SiStripGainFromAsciiFile::getNewObject] detid already exists" << std::endl;
    }
  }
  edm::LogInfo("SiStripGainFromAsciiFile") << ss.str();

  return obj;
}