diff --git a/analyses/pluginMC/MC_JET_IN_HI.cc b/analyses/pluginMC/MC_JET_IN_HI.cc
new file mode 100644
--- /dev/null
+++ b/analyses/pluginMC/MC_JET_IN_HI.cc
@@ -0,0 +1,159 @@
+// -*- C++ -*-
+#include "Rivet/Analysis.hh"
+#include "Rivet/Projections/FinalState.hh"
+#include "Rivet/Projections/FastJets.hh"
+#include "Rivet/Projections/WFinder.hh"
+#include "Rivet/Projections/ZFinder.hh"
+#include "fastjet/tools/Filter.hh"
+#include "fastjet/tools/Pruner.hh"
+#include "Rivet/Tools/AtlasCommon.hh"
+
+namespace Rivet {
+
+  class MC_JET_IN_HI : public Analysis {
+  public:
+
+    /// @name Constructors etc.
+    //@{
+
+    /// Constructor
+    DEFAULT_RIVET_ANALYSIS_CTOR(MC_JET_IN_HI);
+
+    //@}
+
+
+  public:
+
+    string ts(int in) {
+      std::stringstream ss;
+      ss << in;
+      return ss.str();
+    }
+    /// @name Analysis methods
+    //@{
+
+    /// Book histograms and initialise projections before the run
+    void init() {
+      // Declare centrality projection - we use the ATLAS PbPb definition
+      // to be able to compare to data.
+      declareCentrality(ATLAS::SumET_PBPB_Centrality(),"ATLAS_PBPB_CENTRALITY",
+		      "sumETFwd","sumETFwd");
+      // The final state where jets are found
+      FinalState fs(Cuts::abseta < 2.5);
+      declare(fs, "FS");
+ 
+      ZFinder zfinder(fs, Cuts::abseta < 2.5 && Cuts::pT > 30*GeV, PID::MUON, 80*GeV, 100*GeV,
+                      0.2, ZFinder::CLUSTERNODECAY, ZFinder::TRACK);
+      declare(zfinder, "ZFinder");
+
+      // Z+jet jet collections
+      declare(FastJets(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.3), "JetsAK3");
+      declare(FastJets(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.5), "JetsAK5");
+      declare(FastJets(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.7), "JetsAK7");
+      declare(FastJets(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.9), "JetsAK9");
+
+      jetFinders = {"JetsAK3", "JetsAK5", "JetsAK7", "JetsAK9"};
+      
+      for (size_t i = 0; i < jetFinders.size(); ++i) {
+        string s = jetFinders[i];
+	h_zpT.push_back(bookHisto1D(s + "zpT",logspace(50, 1.0,1000)));
+	h_jetpT.push_back(bookHisto1D(s + "jetpT",logspace(50, 1.0,1000)));
+      }
+      incSow = bookCounter("incSow");
+
+      centData = {0., 0.2, 0.4, 0.6, 0.8,};
+      for (size_t i = 0; i < centData.size(); ++i) {
+        c_jetpT[centData[i]] = bookHisto1D("cjetpT" + ts(i),logspace(100, 10.0,1000));
+        c_zpT[centData[i]] = bookHisto1D("czpt" + ts(i),logspace(100, 10.0,1000));
+	sow[centData[i]] = bookCounter("sow_" + ts(i));
+      }
+    }
+
+    bool isBackToBack_zj(const ZFinder& zf, const fastjet::PseudoJet& psjet) {
+      const FourMomentum& z = zf.bosons()[0].momentum();
+      const FourMomentum jmom(psjet.e(), psjet.px(), psjet.py(), psjet.pz());
+      return (deltaPhi(z, jmom) > 7.*M_PI/8. );
+    }
+
+
+    /// Perform the per-event analysis
+    void analyze(const Event& event) {
+      const double weight = event.weight();
+      
+      // Get the Z
+      const ZFinder& zfinder = apply<ZFinder>(event, "ZFinder");
+      if (zfinder.bosons().size() != 1) vetoEvent;
+      Particle z = zfinder.bosons()[0];
+      Particle l1 = zfinder.constituents()[0];
+      Particle l2 = zfinder.constituents()[1];
+
+      // Require a high-pT Z (and constituents)
+      if (l1.pT() < 10*GeV || l2.pT() < 10*GeV || z.pT() < 60*GeV) vetoEvent;
+
+      // Get the centrality
+      const double c = apply<CentralityProjection>(event,"sumETFwd")();
+      auto jetpTItr = c_jetpT.upper_bound(c);
+      auto zpTItr = c_zpT.upper_bound(c);
+      auto sowItr = sow.upper_bound(c);
+      if (jetpTItr == c_jetpT.end() || zpTItr == c_zpT.end() || 
+        sowItr == sow.end()) vetoEvent;
+      sowItr->second->fill(weight);
+      incSow->fill(weight);
+      // Get the  jets
+      for (size_t i = 0; i < jetFinders.size(); ++i ) {
+        const PseudoJets& psjets = apply<FastJets>(event, 
+	  jetFinders[i]).pseudoJetsByPt(30.0*GeV);
+        if (!psjets.empty()) {
+        // Get the leading jet and make sure it's back-to-back with the Z
+        const fastjet::PseudoJet& j0 = psjets[0];
+	if (isBackToBack_zj(zfinder, j0)) {
+	    // Fill the centrality inclusive histograms
+	    h_zpT[i]->fill(z.pT(),weight);
+	    h_jetpT[i]->fill(j0.perp(),weight);
+	    // Fill centrality dept histograms only for R = 0.3
+	    if (i == 0) {
+              jetpTItr->second->fill(j0.perp(),weight);
+	      zpTItr->second->fill(z.pT(),weight);
+	    }
+	  }
+        }
+      }
+    }
+
+
+    /// Normalise histograms etc., after the run
+    void finalize() {
+       for(size_t i = 0; i < jetFinders.size(); ++i) {
+         h_zpT[i]->scaleW(1./incSow->sumW());
+         h_jetpT[i]->scaleW(1./incSow->sumW());
+       }
+       for (size_t i = 0; i < centData.size(); ++i) {
+         c_jetpT[centData[i]]->scaleW(1./sow[centData[i]]->sumW());
+         c_zpT[centData[i]]->scaleW(1./sow[centData[i]]->sumW());
+       }
+    }
+
+
+    //@}
+
+
+  private:
+  vector<string> jetFinders;
+  // Centrality inclusive histograms
+  vector<Histo1DPtr> h_zpT;
+  vector<Histo1DPtr> h_jetpT;
+  CounterPtr incSow;
+  // Centrality intervals
+  vector<double> centData;
+  // Centrality binned histograms
+  map<double, Histo1DPtr> c_jetpT;
+  map<double, Histo1DPtr> c_zpT;
+  map<double, CounterPtr> sow;
+
+  };
+
+
+  // The hook for the plugin system
+  DECLARE_RIVET_PLUGIN(MC_JET_IN_HI);
+
+}
diff --git a/analyses/pluginMC/MC_JET_IN_HI.info b/analyses/pluginMC/MC_JET_IN_HI.info
new file mode 100644
--- /dev/null
+++ b/analyses/pluginMC/MC_JET_IN_HI.info
@@ -0,0 +1,10 @@
+Name: MC_JET_IN_HI
+Summary: Monte Carlo validation, simple jet analysis in heavy ion collisions, using centrality framework.
+Status: VALIDATED
+RunInfo:
+   Heavy Ion event with a Z decaying to muons and an associated jet. Centrality measure taken from ATLAS_2015_PBPBCENTRALITY, which needs to be run, and the corresponding histogram preloaded.
+Authors:
+ - Christian Bierlich <christian.bierlich@thep.lu.se>
+ - Johannes Bellm
+Options:
+ - cent=REF,GEN,IMP,USR