diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt
index 4d432e1..78f81a0 100644
--- a/examples/CMakeLists.txt
+++ b/examples/CMakeLists.txt
@@ -1,26 +1,28 @@
 set(example_dir "${CMAKE_CURRENT_SOURCE_DIR}")
 
 file(GLOB source_files ${example_dir}/*.cc)
 
 if(NOT ROOT_FOUND)
-  list(REMOVE_ITEM source_files "${example_dir}/AnalysisROOT.cc")
+  list(REMOVE_ITEM source_files "${example_dir}/AnalysisROOT.cc" "${example_dir}/WRapidityAnalysisROOT.cc")
 endif()
 
 foreach(source ${source_files})
   get_filename_component(name ${source} NAME_WE)
   add_library(${name}_lib SHARED ${source})
   set_target_properties(${name}_lib PROPERTIES OUTPUT_NAME "${name}")
   if (CMAKE_CXX_COMPILER_ID STREQUAL "Clang")
     # Clang warns about `std::unique_ptr` in C linked (`extern "C"`) function
     # `make_analysis`. However this call is safe from C++, and we don't care
     # about pure C. We can ignore the warning.
     target_compile_options(${name}_lib PRIVATE "-Wno-return-type-c-linkage")
   endif()
   target_compile_options(${name}_lib PRIVATE "-fvisibility=hidden")
   ## link library
   target_link_libraries(${name}_lib HEJ)
 endforeach(source)
 
 if(ROOT_FOUND)
-  target_link_libraries(AnalysisROOT_lib ROOT::Core ROOT::Hist ROOT::RIO)
+  set(ROOT_LIBS ROOT::Core ROOT::Hist ROOT::RIO)
+  target_link_libraries(AnalysisROOT_lib $(ROOT_LIBS))
+  target_link_libraries(WRapidityAnalysisROOT_lib $(ROOT_LIBS))
 endif()
diff --git a/examples/WRapidityAnalysisROOT.cc b/examples/WRapidityAnalysisROOT.cc
new file mode 100644
index 0000000..12e88d3
--- /dev/null
+++ b/examples/WRapidityAnalysisROOT.cc
@@ -0,0 +1,109 @@
+//! Simple HEJ analysis to plot y_W with ROOT
+// To compile analysis:
+// g++ `HEJ-config --cxxflags` `root-config --cflags --libs` `fastjet-config --cxxflags` -fPIC -shared -O2 -fvisibility=hidden -Wl,-soname,libWRapidityAnalysisROOT.so -o libWRapidityAnalysisROOT.so WRapidityAnalysisROOT.cc
+// To run HEJFOG:
+// LD_PRELOAD=/usr/lib64/root/libCore.so.6.26 bin/HEJFOG ../hej/FixedOrderGen/configFO.yml
+#include <memory>     // for std::unique_ptr
+#include <iostream>   // for std::cout
+
+#include "TFile.h"
+#include "TH1.h"
+
+#include "HEJ/Analysis.hh"
+#include "HEJ/Event.hh"
+#include "HEJ/Version.hh"
+
+namespace YAML {
+  class Node;
+}
+
+namespace LHEF {
+  class HEPRUP;
+}
+namespace {
+  void Fill2DHistogram(TH1D & h,double x,double wt)
+  {
+    const int binnumber=h.GetXaxis()->FindBin(x);
+    const double binwidth=h.GetXaxis()->GetBinWidth(binnumber);
+    h.Fill(x,wt/binwidth);
+  }
+
+  class WRapidityAnalysisROOT: public HEJ::Analysis {
+  private:
+    TFile hfile_;
+    TH1D htotal_,hydif_,hyf_,hyb_,hyW_,hyWsampling_;
+  public:
+    WRapidityAnalysisROOT(
+      YAML::Node const & /* config */, LHEF::HEPRUP const & /* heprup */
+    ):hfile_{"HEJ.root","RECREATE",
+             ("Output from "+HEJ::Version::package_name_full()).c_str()},
+      htotal_{"total","total",1,-0.5,.5},
+      hydif_{"ydif","ydif",40,0.,10.},
+      hyf_{"yf","yf",40,-5.,5.},
+      hyb_{"yb","yb",40,-5.,5.},
+      hyW_{"yW","yW",40,-5.,5.},
+      hyWsampling_{"yWsampling","yWsampling",40,-5.,5.}
+    {
+      // Setup histograms
+      htotal_.Sumw2();
+      hydif_.Sumw2();
+      hyf_.Sumw2();
+      hyb_.Sumw2();
+      hyW_.Sumw2();
+    }
+
+    void fill(
+      HEJ::Event const & event,
+      HEJ::Event const & /* FO_event */
+    ) override {
+      const double weight=event.central().weight;
+      htotal_.Fill(0.0,weight);
+      auto const & jets=event.jets(); // jets are sorted in rapidity
+      // find the W momentum
+      auto const & outgoing=event.outgoing();
+      auto the_boson = std::find_if(outgoing.begin(), outgoing.end(), [](auto const & particle)
+                                    { return HEJ::is_AWZ_boson(particle); });
+      Fill2DHistogram(hydif_,
+          fabs(jets.back().rapidity()-jets.front().rapidity()), weight);
+      Fill2DHistogram(hyf_, jets.front().rapidity(), weight);
+      Fill2DHistogram(hyb_, jets.back().rapidity(), weight);
+      Fill2DHistogram(hyW_, the_boson->rapidity(), weight);
+      Fill2DHistogram(hyWsampling_, the_boson->rapidity(),1);
+    }
+
+    bool pass_cuts(
+      HEJ::Event const & event,
+      HEJ::Event const & /* FO_event */
+    ) override {
+      if (event.jets().size()<2)
+        return false;
+      // Find the charged lepton momentum
+      const auto the_decay = event.decays().begin();
+      const auto the_charged_lepton = *std::find_if(
+        the_decay->second.begin(), the_decay->second.end(),
+        [](HEJ::Particle const & p) {
+          return HEJ::is_anylepton(p) && !HEJ::is_anyneutrino(p);
+        }
+      );
+      // momentum: the_charged_lepton.p
+      if (fabs(the_charged_lepton.p.rapidity())>2.5)
+        return false;
+      return true;
+    }
+
+    void finalise() override {
+      hfile_.Write();
+      hfile_.Close();
+      std::cout << "WRapidityAnalysisROOT bids you farewell " << std::endl; // be polite
+    }
+
+  };
+}
+
+extern "C"
+__attribute__((visibility("default")))
+std::unique_ptr<HEJ::Analysis> make_analysis(
+    YAML::Node const & config, LHEF::HEPRUP const & heprup
+){
+  return std::make_unique<WRapidityAnalysisROOT>(config, heprup);
+}