/** \class HLTMuonL2FromL1TPreFilter
 *
 * See header file for documentation
 *
 *  \author J. Alcaraz
 *
 */

#include "HLTMuonL2FromL1TPreFilter.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "DataFormats/Common/interface/Handle.h"
#include "DataFormats/Common/interface/RefToBase.h"
#include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
#include "DataFormats/TrackReco/interface/Track.h"
#include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
#include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
#include "DataFormats/BeamSpot/interface/BeamSpot.h"
#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
#include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
#include "FWCore/Utilities/interface/InputTag.h"

//
// constructors and destructor
//

HLTMuonL2FromL1TPreFilter::HLTMuonL2FromL1TPreFilter(const edm::ParameterSet& iConfig)
    : HLTFilter(iConfig),
      beamSpotTag_(iConfig.getParameter<edm::InputTag>("BeamSpotTag")),
      beamSpotToken_(consumes<reco::BeamSpot>(beamSpotTag_)),
      candTag_(iConfig.getParameter<edm::InputTag>("CandTag")),
      candToken_(consumes<reco::RecoChargedCandidateCollection>(candTag_)),
      previousCandTag_(iConfig.getParameter<edm::InputTag>("PreviousCandTag")),
      previousCandToken_(consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_)),
      seedMapTag_(iConfig.getParameter<edm::InputTag>("SeedMapTag")),
      seedMapToken_(consumes<SeedMap>(seedMapTag_)),
      minN_(iConfig.getParameter<int>("MinN")),
      maxEta_(iConfig.getParameter<double>("MaxEta")),
      absetaBins_(iConfig.getParameter<std::vector<double> >("AbsEtaBins")),
      minNstations_(iConfig.getParameter<std::vector<int> >("MinNstations")),
      minNhits_(iConfig.getParameter<std::vector<int> >("MinNhits")),
      cutOnChambers_(iConfig.getParameter<bool>("CutOnChambers")),
      minNchambers_(iConfig.getParameter<std::vector<int> >("MinNchambers")),
      maxDr_(iConfig.getParameter<double>("MaxDr")),
      minDr_(iConfig.getParameter<double>("MinDr")),
      maxDz_(iConfig.getParameter<double>("MaxDz")),
      min_DxySig_(iConfig.getParameter<double>("MinDxySig")),
      minPt_(iConfig.getParameter<double>("MinPt")),
      nSigmaPt_(iConfig.getParameter<double>("NSigmaPt")),
      matchPreviousCand_(iConfig.getParameter<bool>("MatchToPreviousCand")) {
  using namespace std;

  // check that number of eta bins matches number of nStation cuts
  if (minNstations_.size() != absetaBins_.size() || minNhits_.size() != absetaBins_.size() ||
      (cutOnChambers_ && minNchambers_.size() != absetaBins_.size())) {
    throw cms::Exception("Configuration") << "Number of MinNstations, MinNhits, or MinNchambers cuts "
                                          << "does not match number of eta bins." << endl;
  }

  if (absetaBins_.size() > 1) {
    for (unsigned int i = 0; i < absetaBins_.size() - 1; ++i) {
      if (absetaBins_[i + 1] <= absetaBins_[i])
        throw cms::Exception("Configuration") << "Absolute eta bins must be in increasing order." << endl;
    }
  }

  // dump parameters for debugging
  if (edm::isDebugEnabled()) {
    ostringstream ss;
    ss << "Constructed with parameters:" << endl;
    ss << "    BeamSpotTag = " << beamSpotTag_.encode() << endl;
    ss << "    CandTag = " << candTag_.encode() << endl;
    ss << "    PreviousCandTag = " << previousCandTag_.encode() << endl;
    ss << "    SeedMapTag = " << seedMapTag_.encode() << endl;
    ss << "    MinN = " << minN_ << endl;
    ss << "    MaxEta = " << maxEta_ << endl;
    ss << "    MinNstations = ";
    for (unsigned int j = 0; j < absetaBins_.size(); ++j) {
      ss << minNstations_[j] << " (|eta|<" << absetaBins_[j] << "), ";
    }
    ss << endl;
    ss << "    MinNhits = ";
    for (unsigned int j = 0; j < absetaBins_.size(); ++j) {
      ss << minNhits_[j] << " (|eta|<" << absetaBins_[j] << "), ";
    }
    ss << endl;
    ss << "    CutOnChambers = " << cutOnChambers_ << endl;
    if (cutOnChambers_) {
      ss << "    MinNchambers = ";
      for (unsigned int j = 0; j < absetaBins_.size(); ++j) {
        ss << minNchambers_[j] << " (|eta|<" << absetaBins_[j] << "), ";
      }
    }
    ss << endl;
    ss << "    MaxDr = " << maxDr_ << endl;
    ss << "    MinDr = " << minDr_ << endl;
    ss << "    MaxDz = " << maxDz_ << endl;
    ss << "    MinDxySig = " << min_DxySig_ << endl;
    ss << "    MinPt = " << minPt_ << endl;
    ss << "    NSigmaPt = " << nSigmaPt_ << endl;
    ss << "    saveTags= " << saveTags();
    LogDebug("HLTMuonL2FromL1TPreFilter") << ss.str();
  }
}

HLTMuonL2FromL1TPreFilter::~HLTMuonL2FromL1TPreFilter() = default;

void HLTMuonL2FromL1TPreFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
  edm::ParameterSetDescription desc;
  makeHLTFilterDescription(desc);
  desc.add<edm::InputTag>("BeamSpotTag", edm::InputTag("hltOfflineBeamSpot"));
  desc.add<edm::InputTag>("CandTag", edm::InputTag("hltL2MuonCandidates"));
  desc.add<edm::InputTag>("PreviousCandTag", edm::InputTag(""));
  desc.add<edm::InputTag>("SeedMapTag", edm::InputTag("hltL2Muons"));
  desc.add<int>("MinN", 1);
  desc.add<double>("MaxEta", 2.5);
  desc.add<std::vector<double> >("AbsEtaBins", std::vector<double>(1, 9999.));
  desc.add<std::vector<int> >("MinNstations", std::vector<int>(1, 1));
  desc.add<std::vector<int> >("MinNhits", std::vector<int>(1, 0));
  desc.add<bool>("CutOnChambers", false);
  desc.add<std::vector<int> >("MinNchambers", std::vector<int>(1, 0));
  desc.add<double>("MaxDr", 9999.0);
  desc.add<double>("MinDr", -1.0);
  desc.add<double>("MaxDz", 9999.0);
  desc.add<double>("MinDxySig", -1.0);
  desc.add<double>("MinPt", 0.0);
  desc.add<double>("NSigmaPt", 0.0);
  desc.add<bool>("MatchToPreviousCand", true);
  descriptions.add("hltMuonL2FromL1TPreFilter", desc);
}

//
// member functions
//

// ------------ method called to produce the data  ------------
bool HLTMuonL2FromL1TPreFilter::hltFilter(edm::Event& iEvent,
                                          const edm::EventSetup& iSetup,
                                          trigger::TriggerFilterObjectWithRefs& filterproduct) const {
  // All HLT filters must create and fill an HLT filter object,
  // recording any reconstructed physics objects satisfying (or not)
  // this HLT filter, and place it in the Event.

  using namespace std;
  using namespace edm;
  using namespace reco;
  using namespace trigger;
  //   using namespace l1extra;

  // save Tag
  if (saveTags())
    filterproduct.addCollectionTag(candTag_);

  // get hold of all muon candidates available at this level
  Handle<RecoChargedCandidateCollection> allMuons;
  iEvent.getByToken(candToken_, allMuons);

  // get hold of the beam spot
  Handle<BeamSpot> beamSpotHandle;
  iEvent.getByToken(beamSpotToken_, beamSpotHandle);
  BeamSpot::Point beamSpot = beamSpotHandle->position();

  // get the L2 to L1 map object for this event
  HLTMuonL2ToL1TMap mapL2ToL1(previousCandToken_, seedMapToken_, iEvent);

  // number of eta bins for cut on number of stations
  const std::vector<double>::size_type nAbsetaBins = absetaBins_.size();

  // look at all allMuons,  check cuts and add to filter object
  int n = 0;
  for (auto cand = allMuons->begin(); cand != allMuons->end(); cand++) {
    TrackRef mu = cand->get<TrackRef>();

    // check if this muon passed previous level
    if (matchPreviousCand_ && !mapL2ToL1.isTriggeredByL1(mu))
      continue;

    // eta cut
    if (std::abs(mu->eta()) > maxEta_)
      continue;

    // cut on number of stations
    bool failNstations(false), failNhits(false), failNchambers(false);
    for (unsigned int i = 0; i < nAbsetaBins; ++i) {
      if (std::abs(mu->eta()) < absetaBins_[i]) {
        if (mu->hitPattern().muonStationsWithAnyHits() < minNstations_[i]) {
          failNstations = true;
        }
        if (mu->numberOfValidHits() < minNhits_[i]) {
          failNhits = true;
        }
        if (cutOnChambers_ &&
            (mu->hitPattern().dtStationsWithAnyHits() + mu->hitPattern().cscStationsWithAnyHits() < minNchambers_[i])) {
          failNchambers = true;
        }
        break;
      }
    }
    if (failNstations || failNhits || failNchambers)
      continue;

    //dr cut
    if (std::abs(mu->dxy(beamSpot)) > maxDr_)
      continue;

    //dr cut
    if (std::abs(mu->dxy(beamSpot)) < minDr_)
      continue;

    //dz cut
    if (std::abs(mu->dz(beamSpot)) > maxDz_)
      continue;

    // dxy significance cut (safeguard against bizarre values)
    if (min_DxySig_ > 0 && (mu->dxyError() <= 0 || std::abs(mu->dxy(beamSpot) / mu->dxyError()) < min_DxySig_))
      continue;

    // Pt threshold cut
    double pt = mu->pt();
    double abspar0 = std::abs(mu->parameter(0));
    double ptLx = pt;
    // convert 50% efficiency threshold to 90% efficiency threshold
    if (abspar0 > 0)
      ptLx += nSigmaPt_ * mu->error(0) / abspar0 * pt;
    if (ptLx < minPt_)
      continue;

    // add the good candidate to the filter object
    filterproduct.addObject(
        TriggerMuon, RecoChargedCandidateRef(Ref<RecoChargedCandidateCollection>(allMuons, cand - allMuons->begin())));

    n++;
  }

  // filter decision
  const bool accept(n >= minN_);

  // dump event for debugging
  if (edm::isDebugEnabled()) {
    ostringstream ss;
    ss << "L2mu#" << '\t' << "q*pt" << '\t'  //scientific is too wide
       << '\t' << "q*ptLx" << '\t'           //scientific is too wide
       << '\t' << "eta" << '\t' << "phi" << '\t' << "nStations" << '\t' << "nHits" << '\t' << "dr"
       << '\t'                  //scientific is too wide
       << '\t' << "dz" << '\t'  //scientific is too wide
       << '\t' << "L1seed#" << '\t' << "isPrev" << '\t' << "isFired" << endl;
    ss << "------------------------------------------------------------------------------------------------------------"
          "-----------"
       << endl;
    for (auto cand = allMuons->begin(); cand != allMuons->end(); cand++) {
      TrackRef mu = cand->get<TrackRef>();
      ss << setprecision(2) << cand - allMuons->begin() << '\t' << scientific << mu->charge() * mu->pt() << '\t'
         << scientific
         << mu->charge() * mu->pt() *
                (1. + ((mu->parameter(0) != 0) ? nSigmaPt_ * mu->error(0) / std::abs(mu->parameter(0)) : 0.))
         << '\t' << fixed << mu->eta() << '\t' << fixed << mu->phi() << '\t'
         << mu->hitPattern().muonStationsWithAnyHits() << '\t' << mu->numberOfValidHits() << '\t' << scientific
         << mu->d0() << '\t' << scientific << mu->dz() << '\t' << mapL2ToL1.getL1Keys(mu) << '\t'
         << mapL2ToL1.isTriggeredByL1(mu);
      vector<RecoChargedCandidateRef> firedMuons;
      filterproduct.getObjects(TriggerMuon, firedMuons);
      ss << '\t'
         << (find(firedMuons.begin(),
                  firedMuons.end(),
                  RecoChargedCandidateRef(Ref<RecoChargedCandidateCollection>(allMuons, cand - allMuons->begin()))) !=
             firedMuons.end())
         << endl;
    }
    ss << "------------------------------------------------------------------------------------------------------------"
          "-----------"
       << endl;
    ss << "Decision of filter is " << accept << ", number of muons passing = " << filterproduct.muonSize();
    LogDebug("HLTMuonL2FromL1TPreFilter") << ss.str();
  }

  return accept;
}

// declare this class as a framework plugin
#include "FWCore/Framework/interface/MakerMacros.h"
DEFINE_FWK_MODULE(HLTMuonL2FromL1TPreFilter);