Index: python/DelphesAnalysis/AnalysisEvent.py
===================================================================
--- python/DelphesAnalysis/AnalysisEvent.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/AnalysisEvent.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,277 @@
+from inspect import getargspec
+from ROOT import TLorentzVector,TVector3,TChain,TClass,TDatabasePDG
+from datetime import datetime
+from collections import Iterable
+from types import StringTypes
+from os import path
+
+class AnalysisEvent(TChain):
+   """A class that complements fwlite::Events with analysis facilities.
+      The class provides the following additional functionalities:
+        1. instrumentation for event weight
+             A set of weight classes can be defined, and the event weight
+             is computed and cached using those.
+        2. list of event products used in the analysis
+             It makes the iteration faster by only enabling required branches.
+        3. a list of "producers" of analysis high-level quantities
+             It allows to run "analysis on demand", by automatically running
+             the defined producers to fill the cache, and later use that one.
+        4. a volatile dictionary
+             It allows to use the event as an heterogenous container for
+             any analysis product. The event is properly reset when iterating 
+             to the next event.
+   """
+
+   def __init__(self, inputFiles = '', maxEvents=0):
+     """Initialize the AnalysisEvent like a standard Event, plus additional features."""
+     # initialization of base functionalities
+     TChain.__init__(self,"Delphes","Delphes")
+     if isinstance(inputFiles,Iterable) and not isinstance(inputFiles,StringTypes):
+       for thefile in inputFiles:
+         if path.isfile(thefile): 
+           self.AddFile(thefile)
+         else:
+           print "Warning: ",thefile," do not exist."
+     elif  isinstance(inputFiles,StringTypes):
+       thefile = inputFiles
+       if path.isfile(thefile): 
+         self.AddFile(thefile)
+       else:
+         print "Warning: ",thefile," do not exist."
+     else:
+       print "Warning: invalid inputFiles"
+     self.BuildIndex("Event[0].Number")
+     self.SetBranchStatus("*",0)
+     self._eventCounts        = 0
+     self._maxEvents          = maxEvents
+     # additional features:
+     # 1. instrumentation for event weight
+     self._weightCache = {}
+     self._weightEngines = {}
+     # 2. a list of event products used in the analysis
+     self._collections = {} 
+     self._branches = dict((b,False) for b in map(lambda b:b.GetName(),self.GetListOfBranches()))
+     # 3. a list of "producers" of analysis high-level quantities
+     self._producers = {}
+     # 4. volatile dictionary. User can add any quantity to the event and it will be 
+     #    properly erased in the iteration step. 
+     self.__dict__["vardict"] = {}
+     
+   def addWeight(self, name, weightClass):
+     """Declare a new class (engine) to compute the weights.
+        weightClass must have a weight() method returning a float."""
+     if name in self._weightEngines:
+       raise KeyError("%s weight engine is already declared" % name)
+     self._weightEngines[name] = weightClass
+     self._weightCache.clear()
+
+   def delWeight(self, name):
+     """Remove one weight engine from the internal list."""
+     # just to clean the dictionnary
+     del self._weightEngines[name]
+     self._weightCache.clear()
+
+   def weight(self, weightList=None, **kwargs):
+     """Return the event weight. Arguments:
+         * weightList is the list of engines to use, as a list of strings. 
+              Default: all defined engines.
+         * the other named arguments are forwarded to the engines.
+        The output is the product of the selected individual weights."""
+     # first check in the cache if the result is there already
+     if weightList is None:
+       weightList=self._weightEngines.keys()
+     kwargs["weightList"] = weightList
+     # compute the weight or use the cached value
+     myhash = self._dicthash(kwargs)
+     if not myhash in self._weightCache : 
+       w = 1.
+       for weightElement in weightList:
+         engine = self._weightEngines[weightElement]
+         engineArgs = getargspec(engine.weight).args
+         subargs = dict((k,v) for k,v in kwargs.iteritems() if k in engineArgs)
+         w *= self._weightCache.setdefault("weightElement:%s # %s" %(weightElement,self._dicthash(subargs)),engine.weight(self,**subargs))
+       self._weightCache[myhash] = w
+     return self._weightCache[myhash]
+
+   def addCollection(self, name, inputTag):
+     """Register an event collection as used by the analysis.
+        Example: addCollection("myjets","jets") 
+        Note that the direct access to the branch is still possible but unsafe."""
+     if name in self._collections:
+       raise KeyError("%r collection is already declared", name)
+     if name in self._producers:
+       raise KeyError("%r is already declared as a producer", name)
+     if hasattr(self,name):
+       raise AttributeError("%r object already has attribute %r" % (type(self).__name__, name))
+     if inputTag not in self._branches:
+       raise AttributeError("%r object has no branch %r" % (type(self).__name__, inputTag))
+     self._collections[name] = inputTag
+     self.SetBranchStatus(inputTag+"*",1)
+     self._branches[inputTag] = True
+
+   def removeCollection(self,name):
+     """Forget about the named event collection.
+        This method will delete both the product from the cache (if any) and the definition.
+        To simply clear the cache, use "del event.name" instead. """
+     self.SetBranchStatus(self._collections[name]+"*",0)
+     self._branches[self._collections[name]] = False
+     del self._collections[name]
+     if name in self.vardict:
+       delattr(self,name)
+
+   def getCollection(self,name):
+     """Retrieve the event product or return the cached collection.
+        Note that the prefered way to get the collection is instead to access the "event.name" attribute."""
+     if not name in self._collections:
+       raise AttributeError("%r object has no attribute %r" % (type(self).__name__,name))
+     if not name in self.vardict:
+       self.vardict[name] = TChain.__getattr__(self,self._collections[name])
+     return getattr(self,name)
+
+   def addProducer(self,name,producer,**kwargs):
+     """Register a producer to create new high-level analysis objects."""
+     # sanity checks
+     if name in self._producers:
+       raise KeyError("%r producer is already declared", name)
+     if name in self._collections:
+       raise KeyError("%r is already declared as a collection", name)
+     if hasattr(self,name):
+       raise AttributeError("%r object already has attribute %r" % (type(self).__name__, attr))
+     # remove name and producer from kwargs
+     if "name" in kwargs: del kwargs["name"]
+     if "producer" in kwargs: del kwargs["producer"]
+     # store
+     self._producers[name] = (producer,kwargs)
+
+   def removeProducer(self,name):
+     """Forget about the producer.
+        This method will delete both the product from the cache (if any) and the producer.
+        To simply clear the cache, use "del event.name" instead."""
+     del self._producers[name]
+     if name in self.vardict:
+       delattr(self,name)
+
+   def event(self):
+     """Event number"""
+     if self._branches["Event"]:
+       return self.Event.At(0).Number
+     else:
+       return 0
+
+   def to(self,event):
+     """Jump to some event"""
+     self.GetEntryWithIndex(event)
+ 
+   def __getitem__(self,index):
+     """Jump to some event"""
+     self.GetEntryWithIndex(index)
+     return self
+
+   def __iter__ (self):
+     """Iterator"""
+     self._eventCounts = 0
+     while self.GetEntry(self._eventCounts):
+       self.vardict.clear()
+       self._weightCache.clear()
+       yield self
+       self._eventCounts += 1
+       if self._maxEvents > 0 and self._eventCounts >= self._maxEvents:
+         break
+
+   def __getattr__(self, attr):
+     """Overloaded getter to handle properly:
+          - volatile analysis objects
+          - event collections
+          - data producers"""
+     if attr in self.__dict__["vardict"]:
+       return self.vardict[attr]
+     if attr in self._collections:
+       return self.vardict.setdefault(attr, TChain.__getattr__(self,self._collections[attr]))
+     if attr in self._producers:
+       return self.vardict.setdefault(attr, self._producers[attr][0](self, **self._producers[attr][1]))
+     return TChain.__getattr__(self,attr)
+
+   def __setattr__(self, name, value):
+     """Overloaded setter that puts any new attribute in the volatile dict."""
+     if name in self.__dict__  or not "vardict" in self.__dict__ or name[0]=='_':
+       self.__dict__[name] = value
+     else:
+       if name in self._collections or name in self._producers:
+         raise AttributeError("%r object %r attribute is read-only (event collection)" % (type(self).__name__, name))
+       self.vardict[name] = value
+
+   def __delattr__(self, name):
+     """Overloaded del method to handle the volatile internal dictionary."""
+     if name=="vardict":
+       raise AttributeError("%r object has no attribute %r" % (type(self).__name__, name))
+     if name in self.__dict__:
+       del self.__dict__[name]
+     elif name in self.vardict:
+       del self.vardict[name]
+     else:
+       raise AttributeError("%r object has no attribute %r" % (type(self).__name__, name))
+
+   def _dicthash(self,dict):
+     return (lambda d,j='=',s=';': s.join([j.join((str(k),str(v))) for k,v in d.iteritems()]))(dict)
+
+   def __str__(self):
+     """Event text dump."""
+     dictjoin = lambda d,j=' => ',s='\n': s.join([j.join((str(k),str(v))) for k,v in d.iteritems()])
+     mystring = "=================================================================\n"
+     # general information
+     if self._branches["Event"]:
+       mystring += str(self.Event.At(0))
+     else:
+       mystring += "Event %d\n" % self.GetReadEvent()
+     mystring += "-----------------------------------------------------------------\n"
+     # weights
+     if len(self._weightCache)==0: 
+       mystring += "No weight computed so far. Default weight is %f.\n" % self.weight()
+     else:
+       mystring += "Weights:\n"
+       mystring += dictjoin(self._weightCache)
+     mystring += "\n-----------------------------------------------------------------\n"
+     # list the collections
+     mystring += "Collections:\n"
+     for colname in self._collections.keys():
+       collection = self.getCollection(colname)
+       if collection.GetEntries()>0:
+         if collection.At(0).IsA()==TClass.GetClass("HepMCEvent"):
+           pass
+         else:
+           mystring += "*** %s has %d element(s)\n" % (colname,collection.GetEntries())
+           mystring += reduce(lambda a,b: a+b,map(str,collection))
+     mystring += "\n-----------------------------------------------------------------\n"
+     # list the registered producers
+     mystring += "Producers:\n"
+     mystring += dictjoin(self._producers)
+     mystring += "\n-----------------------------------------------------------------\n"
+     # list the content of vardict, excluding collections
+     mystring += "Content of the cache:\n"
+     for k, v in self.vardict.iteritems():
+       if k in self._collections.keys() : continue
+       if isinstance(v,Iterable) and not isinstance(v,StringTypes):
+         try:
+           thisstring = "%s => vector of %d objects(s)\n" % (k,len(v))
+         except:
+           mystring += "%s => %s\n"%(k,str(v))
+         else:
+           try:
+             for it,vec in enumerate(v):
+               thisstring += "%s[%d] = %s\n"%(k,it,str(vec))
+           except: 
+             mystring += "%s => %s\n"%(k,str(v))
+           else:
+             mystring += thisstring
+       else:
+         mystring += "%s => %s\n"%(k,str(v))
+     return mystring
+
+   def decayTree(self, genparticles):
+     db = TDatabasePDG()
+     theString = ""
+     for part in genparticles:
+       if part.M1==-1 and part.M2==-1:
+         theString += part.printDecay(db, genparticles)
+     return theString
+
Index: python/DelphesAnalysis/BaseControlPlots.py
===================================================================
--- python/DelphesAnalysis/BaseControlPlots.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/BaseControlPlots.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,166 @@
+#! /usr/bin/env python
+import ROOT 
+from CPconfig import configuration
+
+def getArgSet(controlplots):
+  assert isinstance(controlplots,list)
+  merged = ROOT.RooArgSet()
+  for ctrlplot in controlplots:
+    merged.add(ctrlplot._obsSet)
+  return merged
+
+def runTest(path='../testfiles/', controlPlots=None):
+  # these modules are only needed in that function, used for debugging.
+  # so we only import them here.
+  import Delphes
+  import os
+  import AnalysisEvent
+  import EventSelection
+  assert isinstance(controlPlots, BaseControlPlots)
+  if os.path.isdir(path):
+    dirList=os.listdir(path)
+    files=[]
+    for fname in dirList:
+      files.append(path+fname)
+  elif os.path.isfile(path):
+    files=[path]
+  else:
+    files=[]
+  events = AnalysisEvent.AnalysisEvent(files)
+  EventSelection.prepareAnalysisEvent(events)
+  controlPlots.beginJob()
+  i = 0
+  for event in events:
+    if i%100==0 : print "Processing... event ", i
+    controlPlots.processEvent(event)
+    i += 1
+  controlPlots.endJob()
+
+class BaseControlPlots:
+    """A class to create control plots"""
+    
+    def __init__(self, dir=None, purpose="generic", dataset=None, mode="plots"):
+      """Initialize the ControlPlots, creating output file if needed. If no file is given, it means it is delegated."""
+      self._mode = mode
+      self._purpose = purpose
+      # for plots
+      if self._mode=="plots":
+        # create output file if needed. If no file is given, it means it is delegated
+        if dir is None:
+          self._f = ROOT.TFile(configuration.defaultFilename+".root", "RECREATE")
+          self._dir = self._f.mkdir(purpose)
+        else :
+          self._f = None
+          self._dir = dir
+        self._h_vector = { }
+      # for ntuples
+      if self._mode=="dataset":
+        self._obsSet = ROOT.RooArgSet()
+        if dataset is None:
+          self._rds = ROOT.RooDataSet(configuration.RDSname, configuration.RDSname, self._obsSet)
+	  self._ownedRDS = True
+	else:
+	  self._rds = dataset
+	  self._ownedRDS = False
+        self._rrv_vector = { }
+        self._rooCategories = { }
+    
+    def beginJob(self):
+      """Declare histograms, and for derived classes instantiate handles. Must be overloaded.""" 
+      raise NotImplementedError
+
+    def defineCategories(self, categories):
+      """Define the categories, given a list of names. Only works for datasets"""
+      if self._mode!="dataset": return
+      for i, name in enumerate(categories):
+        rc = ROOT.RooCategory("rc_"+self._purpose+"_"+str(i),name.split('/')[-1])
+	rc.defineType("not_acc",0)
+	rc.defineType("acc",1)
+	self._rooCategories[i] = rc
+	self._rds.addColumn(rc)
+	self._obsSet.add(rc)
+
+    def addHisto(self,*args):
+      """Add one histograms to the list of products. Arguments are as for TH1F."""
+      # this fills a distionnary name <-> histogram
+      self._dir.cd()
+      self._h_vector[args[0]] = ROOT.TH1F(*args)
+
+    def addVariable(self,*args):
+      """Add one variable to the list of products. Arguments are as for RooRealVar."""
+      # this fills a distionnary name <-> RooRealVar
+      self._rrv_vector[args[0]] = ROOT.RooRealVar(self._purpose+args[0],*(args[1:]))
+      self._obsSet.add(self._rrv_vector[args[0]])
+      self._rds.addColumn(self._rrv_vector[args[0]])
+    
+    def add(self, *args):
+      """Add one item to the list of products. Arguments are as for TH1F."""
+      if self._mode=="plots":
+        self.addHisto(*args)
+      else:
+        self.addVariable(*[args[i] for i in [0,1,3,4]])
+
+    def processEvent(self, event, weight = 1.):
+      """process event and fill histograms"""
+      self.fill(self.process(event),weight)
+     
+    def process(self, event):
+      """Process event data to extract histogrammables. Must be overloaded."""
+      raise NotImplementedError
+      # this is just an example, we never get there
+      # that method must return a dictionnary name <-> value that matches self._h_vector
+      result = { }
+      result["var1"] = 1.
+      result["var2"] = 2.3
+      result["var3"] = 5.711
+      return result
+
+    def setCategories(self, categories):
+      """Set the categories, given a list of booleans. Only works for datasets"""
+      if self._mode!="dataset": return
+      for c, flag in enumerate(categories):
+        if flag: self._rooCategories[c].setIndex(1)
+	else: self._rooCategories[c].setIndex(0)
+
+    def fillPlots(self, data, weight = 1.):
+      """Fills histograms with the data provided as input."""
+      for name,value in data.items():
+        if isinstance(value,list):
+          for val in value: self._h_vector[name].Fill(val,weight)
+        else:
+          self._h_vector[name].Fill(value,weight)
+
+    def fillRDS(self, data):
+      """Fills roodataset with the data provided as input."""
+      for rrv in self._rrv_vector:
+        # data is not guaranteed to contain all variables,
+        # so if we don't do this, these variables will keep the old value
+        self._rrv_vector[rrv].setVal(-1000)
+      for name,value in data.items():
+        if isinstance(value,list):
+	  #for now, we only store scalars, not vectors
+	  pass
+        else:
+	  self._rrv_vector[name].setVal(value)
+      if self._ownedRDS:
+        self._rds.add(self._obsSet)  
+
+    def fill(self, data, weight = 1.):
+      """Fills whatever must be filled in"""
+      if self._mode=="plots":
+        self.fillPlots(data, weight)
+      else:
+        self.fillRDS(data)
+
+    def endJob(self):
+      """Save and close."""
+      if self._mode=="plots":
+        self._dir.cd()
+        self._dir.Write()
+        if not self._f is None:
+          self._f.Close()
+      else:
+        if self._ownedRDS:
+          ws  = ROOT.RooWorkspace(self._purpose,self._purpose)
+          getattr(ws,'import')(self._rds) 
+          ws.writeToFile(configuration.defaultFilename+"_"+self._purpose+".root") 
Index: python/DelphesAnalysis/BaseWeightClass.py
===================================================================
--- python/DelphesAnalysis/BaseWeightClass.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/BaseWeightClass.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,10 @@
+# Every weight class must contain a weight method, taking as first argument the event.
+# Additional methods and arguments are allowed.
+# Note that weight classes do not have to derive from BaseWeightClass.
+
+class BaseWeightClass:
+  """A class to reweight AnalysisEvents"""
+  def weight( self, event):
+     """Lepton eff weight"""
+     raise NotImplementedError
+
Index: python/DelphesAnalysis/CPconfig.py
===================================================================
--- python/DelphesAnalysis/CPconfig.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/CPconfig.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,41 @@
+#configuration of the ControlPlot machinery
+
+from collections import namedtuple
+controlPlot     = namedtuple("controlPlot",    ["label","module","classname","kwargs"])
+eventCollection = namedtuple("eventCollection",["label","collection"])
+eventProducer   = namedtuple("eventProducer",  ["label","module","function","kwargs"])
+eventWeight     = namedtuple("eventWeight",    ["label","module","classname","kwargs"])
+
+class configuration:
+  # default I/O
+  defaultFilename = "controlPlots"
+  RDSname = "rds_delphes"
+  WSname = "workspace_ras"
+
+  # mode: plots or dataset
+  runningMode = "plots"
+
+  # event selection class
+  eventSelection = ""
+
+  # control plot classes
+  controlPlots = [ ]
+
+  # event content: lists of eventCollection, eventProducer, and eventWeight objects respectively.
+  eventCollections = [ ]
+  eventProducers   = [ ]
+  eventWeights     = [ ]
+
+class eventDumpConfig:
+  # fine-tuning of the event content for display
+  productsToPrint   = [ ] # list of product to display (use the producer label)
+  collectionsToHide = [ ] # collections used in the analysis but not printed (use the collection label) 
+
+# import the actual implementation of the configuration
+import os
+theConfig = os.getenv("DelphesAnalysisCfg")
+if theConfig is not None:
+  configImplementation = __import__(os.path.splitext(theConfig)[0])
+  configuration = configImplementation.configuration
+  eventDumpConfig = configImplementation.eventDumpConfig
+
Index: python/DelphesAnalysis/ControlPlots.py
===================================================================
--- python/DelphesAnalysis/ControlPlots.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/ControlPlots.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,211 @@
+#!/usr/bin/env python 
+
+#######################################################################################
+###  Option parsing and main routine  #################################################
+#######################################################################################
+
+from optparse import OptionParser
+import sys
+import os
+usage="""%prog [options]"""
+description="""A simple script to generate control plots."""
+epilog="""Example:
+./ControlPlots.py -i ./ -o controlPlots_demo.root --all
+"""
+parser = OptionParser(usage=usage,add_help_option=True,description=description,epilog=epilog)
+parser.add_option("-c", "--conf", dest="conf", default=None,
+                  help="Configuration file.")
+parser.add_option("-i", "--inputPath", dest="path",
+                  help="Read input file from DIR.", metavar="DIR")
+parser.add_option("-o", "--output", dest='outputname', default=None,
+                  help="Save output as FILE.", metavar="FILE")
+parser.add_option("--all",action="store_true",dest="all",
+                  help="Process all levels.")
+parser.add_option("-l", "--level", dest="levels",
+                  help="Specify a coma-separated list of levels to be processed. No space is allowed.")
+parser.add_option("--Njobs", type="int", dest='Njobs', default="1",
+                  help="Number of jobs when splitting the processing.")
+parser.add_option("--jobNumber", type="int", dest='jobNumber', default="0",
+                  help="Number of the job is a splitted set of jobs.")
+
+(options, args) = parser.parse_args()
+
+#special treatment of the config file... the rest of the options will be parsed in main.
+if options.conf is not None:
+  try:
+    theUserConf = __import__(os.path.splitext(options.conf)[0])
+  except:
+    raise ImportError("%s is not a valid configuration file."%options.conf)
+  else:
+    os.environ["DelphesAnalysisCfg"]=options.conf
+from CPconfig import configuration
+if options.outputname is None:
+  options.outputname = configuration.defaultFilename+".root"
+
+import Delphes
+import ROOT
+import itertools
+import time
+from importlib import import_module
+from AnalysisEvent import AnalysisEvent
+from BaseControlPlots import getArgSet
+import EventSelection
+import cProfile
+
+def main(options):
+  """simplistic program main"""
+  # do basic arg checking
+  if options.path is None: 
+    print "Error: no input path specified."
+    parser.print_help()
+    return
+  levels = []
+  if options.all:
+    levels = range(EventSelection.eventCategories())
+  elif not options.levels is None:
+    levels= map(int,options.levels.split(','))
+  levels.sort()
+  if len(levels)==0:
+    print "Error: no level specified for processing."
+    parser.print_help()
+    return
+  if min(levels)<0:
+    print "Error: levels must be positive integers."
+    parser.print_help()
+    return
+  if max(levels)>=EventSelection.eventCategories():
+    print "Error: last level is",EventSelection.eventCategories()-1
+    parser.print_help()
+    return
+  if options.Njobs<1:
+    print "Error: Njobs must be strictly positive."
+    parser.print_help()
+    return
+  if options.jobNumber>=options.Njobs:
+    print "Error: jobNumber must be strictly smaller than Njobs."
+    parser.print_help()
+    return
+  # if all ok, run the procedure
+  runAnalysis(path=options.path,outputname=options.outputname, levels=levels, Njobs=options.Njobs, jobNumber=options.jobNumber)
+
+#######################################################################################
+### Central Routine: manage input/output, loop on events, manage weights and plots  ###
+#######################################################################################
+
+def runAnalysis(path, levels, outputname="controlPlots.root", Njobs=1, jobNumber=1):
+  """produce all the plots in one go"""
+
+  # inputs
+  if os.path.isdir(path):
+    dirList=list(itertools.islice(os.listdir(path), jobNumber, None, Njobs))
+    files=[]
+    for fname in dirList:
+      files.append(path+"/"+fname)
+  elif os.path.isfile(path):
+    files=[path]
+  else:
+    files=[]
+
+  # output
+  output = ROOT.TFile(outputname, "RECREATE")
+  if configuration.runningMode=="dataset":
+    ROOT.RooAbsData.setDefaultStorageType(ROOT.RooAbsData.Tree)
+    rds = ROOT.RooDataSet(configuration.RDSname, configuration.RDSname, ROOT.RooArgSet())
+
+  # events iterator, plus configuration of standard collections and producers
+  events = AnalysisEvent(files)
+  EventSelection.prepareAnalysisEvent(events)
+
+  # prepare the plots
+  controlPlots=[]
+  if configuration.runningMode=="plots":
+    leafList = [None]*EventSelection.eventCategories()
+    createDirectory(EventSelection.categoriesHierarchy(), output, leafList)
+    for levelDir in leafList:
+      levelPlots=[]
+      for cp in configuration.controlPlots:
+        levelPlots.append(getattr(import_module(cp.module),cp.classname)(dir=levelDir.mkdir(cp.label),mode="plots"))
+      controlPlots.append(levelPlots)
+  else:
+    for cp in configuration.controlPlots:
+      controlPlots.append(getattr(import_module(cp.module),cp.classname)(dir=None, mode="dataset", dataset=rds))
+
+  # book histograms (separate iteration for clarity)
+  if configuration.runningMode=="plots":
+    for level in levels:
+     for conf,cp in zip(configuration.controlPlots,controlPlots[level]):
+       cp.beginJob(**conf.kwargs)
+  else:
+    for conf,cp in zip(configuration.controlPlots,controlPlots):
+      cp.beginJob(**conf.kwargs)
+    for cp in controlPlots[:1]:
+      cp.defineCategories(EventSelection.categoryNames)
+
+  # process events
+  i = 0
+  t0 = time.time()
+  for event in events:
+    # printout
+    if i%100==0 : 
+      print "Processing... event %d. Last batch in %f s." % (i,(time.time()-t0))
+      t0 = time.time()
+    if configuration.runningMode=="plots":
+      # loop on channels
+      plots = filter(lambda x: EventSelection.isInCategory(x,event.category) ,levels)
+      # process the event once (for the first level)
+      selectionPlotsData=[]
+      for level in plots[:1]:
+        for cp in controlPlots[level]:
+          selectionPlotsData.append(cp.process(event))
+      # fill the histograms
+      for level in plots:
+        for cp, data in zip(controlPlots[level],selectionPlotsData):
+          cp.fill(data, event.weight(category=level))
+    else:
+      for cp in controlPlots[:1]:
+        # set categories (first CP only)
+        cp.setCategories(map(lambda c:EventSelection.isInCategory(c, event.category),range(EventSelection.eventCategories())))
+      for cp in controlPlots:
+        # process event (all CP)
+        cp.processEvent(event)
+      # add to the dataset
+      rds.add(getArgSet(controlPlots))
+    i += 1
+
+  # save all
+  if configuration.runningMode=="plots":
+    for level in levels:
+     for cp in controlPlots[level]: 
+       cp.endJob()
+  else:
+   for cp in controlPlots: 
+     cp.endJob()
+
+  # for dataset, write the merged RDS to file
+  if configuration.runningMode=="dataset":
+    output.cd()
+    ws_ras = ROOT.RooWorkspace(configuration.WSname, configuration.WSname)
+    getattr(ws_ras,'import')(rds.get())
+    output.Add(ws_ras)
+    ws_ras.Write()
+    rds.tree().Write()
+  
+  # close the file
+  output.Close()
+
+def createDirectory(dirStructure, directory, leafList):
+  """Recursively creates the directories for the various stages"""
+  for key,item in dirStructure.iteritems():
+    if isinstance(item, dict):
+      createDirectory(item, directory.mkdir(key), leafList)
+    else:
+      leafList[item]=directory.mkdir("stage_"+str(item),key)
+
+#######################################################################################
+### Program bootstrap  ################################################################
+#######################################################################################
+
+if __name__ == "__main__":
+  main(options)
+  #cProfile.run('main(options)', 'controlPlots.prof')
+
Index: python/DelphesAnalysis/Delphes.py
===================================================================
--- python/DelphesAnalysis/Delphes.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/Delphes.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,160 @@
+import ROOT
+import libPyROOT as _root
+ROOT.gSystem.Load("libDelphes")
+
+#####################################################
+### Definition of the string conversion methods   ###
+#####################################################
+
+def _Event__str__( self ):
+  theString = "Event number %i\n" % self.Number
+  theString += "Read in %f and processed in %f.\n" % (self.ReadTime, self.ProcTime)
+  return theString
+
+_root.MakeRootClass( "Event" ).__str__    = _Event__str__
+
+def _LHCOEvent__str__( self ):
+  theString = "Trigger word: %b\n" % self.Trigger
+  return theString
+
+_root.MakeRootClass( "LHCOEvent" ).__str__    = _LHCOEvent__str__
+
+def _LHEFEvent__str__( self ):
+  hepmcstr = "Process ID: %i\n" % self.ProcessID
+  hepmcstr += "ScalePDF = %f\n" % self.ScalePDF
+  hepmcstr += "Alpha QED = %f\n" % self.AlphaQED
+  hepmcstr += "Alpha DCQ = %f\n" % self.AlphaQCD
+  hepmcstr += "Event weight = %f\n" % self.Weight
+  return hepmcstr
+
+_root.MakeRootClass( "LHEFEvent" ).__str__    = _LHEFEvent__str__
+
+def _HepMCEvent__str__( self ):
+  hepmcstr = "Process ID: %i\n" % self.ProcessID
+  hepmcstr += "Energy scale: %f\n" % self.Scale
+  hepmcstr += "Alpha QED = %f\n" % self.AlphaQED
+  hepmcstr += "Alpha QCD = %f\n" % self.AlphaQCD
+  hepmcstr += "Number of multi parton interactions: %f\n" % self.MPI
+  hepmcstr += """Q-scale used in evaluation of PDFs (in GeV): %f\n""" % self.ScalePDF
+  hepmcstr += """Fraction of beam momentum carried by %i parton ("beam side"): %f\n""" % (self.ID1,self.X1)
+  hepmcstr += """Fraction of beam momentum carried by %i parton ("target side"): %f\n""" % (self.ID2,self.X2)
+  hepmcstr += "PDF (id1, x1, Q) = %f\n" % self.PDF1 
+  hepmcstr += "PDF (id2, x2, Q) = %f\n" % self.PDF2 
+  return hepmcstr
+
+_root.MakeRootClass( "HepMCEvent" ).__str__    = _HepMCEvent__str__
+
+def _GenParticle__str__( self ):
+  thestring =  "Particle HEP ID number = %d\n" % self.PID
+  thestring += "Particle status: %d\n" % self.Status
+  if self.IsPU: thestring += "That particle is from PU.\n"
+  thestring += "Particle 1st mother: %d\n" % self.M1
+  thestring += "Particle 2nd mother: %d\n" % self.M2
+  thestring += "Particle 1st daughter: %d\n" % self.D1
+  thestring += "Particle last daughter: %d\n" % self.D2
+  thestring += "(Px,Py,Pz,E) = (%f,%f,%f,%f)\n" % (self.Px,self.Py,self.Pz,self.E)
+  thestring += "(Pt,Eta,Phi) = (%f,%f,%f)\n" % (self.PT,self.Eta,self.Phi)
+  thestring += "Rapidity = %f\n" % self.Rapidity
+  thestring += "Charge: %d\n" % self.Charge
+  thestring += "Mass: %f\n" % self.Mass
+  thestring += "Vertex position: (x,y,z,t) = (%f,%f,%f,%f)\n" % (self.X,self.Y,self.Z,self.T)
+  return thestring
+
+_root.MakeRootClass( "GenParticle" ).__str__    = _GenParticle__str__
+
+def _GenParticle_printDecay( self, db, particles, pre="" ):
+  thestring = "%s (%f,%f,%f,%f) (status %d)\n" % (db.GetParticle(self.PID).GetName(),self.Px,self.Py,self.Pz,self.E, self.Status)
+  for daughter in range(self.D1,self.D2+1):
+    if daughter==-1: break
+    if daughter<self.D2 : thestring += "%s+->%s" % (pre,particles.At(daughter).printDecay(db,particles,pre+"|   "))
+    else: thestring += "%s+->%s" % (pre,particles.At(daughter).printDecay(db,particles,pre+"    "))
+  return thestring
+
+_root.MakeRootClass( "GenParticle" ).printDecay = _GenParticle_printDecay
+
+def _MissingET__str__( self ):
+  thestring =  "Mising transverse energy: %f\n" % self.MET
+  thestring += "Mising energy azimuthal angle: %f\n" % self.Phi
+  return thestring
+ 
+_root.MakeRootClass( "MissingET" ).__str__    = _MissingET__str__
+
+def _ScalarHT__str__( self ):
+  thestring = "Scalar sum of transverse momenta: %f\n" %self.HT
+
+_root.MakeRootClass( "ScalarHT" ).__str__    = _ScalarHT__str__
+
+def _Photon__str__( self ):
+  thestring = "(Pt,Eta,Phi) =(%f,%f,%f)\n" % (self.PT,self.Eta,self.Phi)
+  return thestring
+
+_root.MakeRootClass( "Photon" ).__str__    = _Photon__str__
+
+def _Electron__str__( self ):
+  thestring = "(Pt,Eta,Phi) =(%f,%f,%f)\n" % (self.PT,self.Eta,self.Phi)
+  thestring += "charge: %d\n" % self.Charge
+  return thestring 
+
+_root.MakeRootClass( "Electron" ).__str__    = _Electron__str__
+
+def _Muon__str__( self ):
+  thestring = "(Pt,Eta,Phi) =(%f,%f,%f)\n" % (self.PT,self.Eta,self.Phi)
+  thestring += "charge: %d\n" % self.Charge
+  return thestring 
+
+_root.MakeRootClass( "Muon" ).__str__    = _Muon__str__
+
+def _Jet__str__( self ):
+  thestring = "(Pt,Eta,Phi) =(%f,%f,%f)\n" % (self.PT,self.Eta,self.Phi)
+  thestring += "Charge: %d\n" % self.Charge
+  thestring += "Mass: %f\n" % self.Mass
+  if self.BTag: thestring += "This jet is b-tagged.\n"
+  if self.TauTag: thestring += "This jet is tagged as a tau-jet.\n" 
+  thestring += "Jet radius in (eta,phi) = (%f,%f)\n" % (self.DeltaEta,self.DeltaPhi)
+  thestring += "Number of constituents: %d\n" % self.Constituents.GetEntries()
+  thestring += "H/E ratio in the calorimeters: %f\n" % self.EhadOverEem
+  nch = 0
+  ech = 0.
+  for constit in self.Constituents:
+    if constit.Charge != 0:
+      fourvector = ROOT.TLorentzVector()
+      fourvector.SetPtEtaPhiM(constit.PT,constit.Eta,constit.Phi,0.)
+      nch += 1
+      ech += fourvector.E()
+  thestring += "Charged Energy fraction: %f\n" % (ech/self.P4().E())
+  thestring += "Charged Multiplicity: %d\n" % nch
+  return thestring
+
+_root.MakeRootClass( "Jet" ).__str__    = _Jet__str__
+
+def _Track__str__( self ):
+  thestring =  "Particle HEP ID number = %d\n" % self.PID
+  thestring += "(Pt,Eta,Phi) = (%f,%f,%f)\n" % (self.PT,self.Eta,self.Phi)
+  thestring += "Charge: %d\n" % self.Charge
+  thestring += "At calorimeter surface: (Pt,Eta,Phi) = (%f,%f,%f)\n" % (self.PT,self.EtaOuter,self.PhiOuter)
+  thestring += "Vertex position: (x,y,z) = (%f,%f,%f)\n" % (self.X,self.Y,self.Z)
+  thestring += "Position at calorimer surface: (x,y,z) = (%f,%f,%f)\n" % (self.XOuter,self.YOuter,self.ZOuter)
+  return thestring 
+
+_root.MakeRootClass( "Track" ).__str__    = _Track__str__
+
+def _Tower__str__( self):
+  thestring  = "(Et,Eta,Phi) = (%f,%f,%f)\n" % (self.ET,self.Eta,self.Phi)
+  thestring += "(E,Eem,Ehad) = (%f,%f,%f)\n" % (tower.E,tower.Eem,tower.Ehad)
+  thestring += "Edges: %f, %f, %f, %f\n" % (tower.Edges[0],tower.Edges[1],tower.Edges[2],tower.Edges[3])
+  return thestring
+
+_root.MakeRootClass( "Tower" ).__str__    = _Tower__str__
+
+#####################################################
+### Definition of additional ROOT methods         ###
+#####################################################
+
+def _lorentzVector__str__( self ):
+  theString =  "(pt, eta, phi) = (%f,%f,%f)\n" % (self.Pt(), self.Eta(), self.Phi())
+  theString += "mass = %f, p = %f, mt = %f\n"  % (self.M(), self.P(), self.Mt())
+  return theString
+
+_root.MakeRootClass( "TLorentzVector" ).__str__    = _lorentzVector__str__
+
+
Index: python/DelphesAnalysis/DumpEventInfo.py
===================================================================
--- python/DelphesAnalysis/DumpEventInfo.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/DumpEventInfo.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,39 @@
+#!/usr/bin/env python 
+import os
+import ROOT
+import Delphes
+from AnalysisEvent import AnalysisEvent
+import EventSelection
+from CPconfig import configuration, eventDumpConfig
+
+def DumpEventInfo(event=None, eventNumber=None, path="./"):
+  """Dump informations about a given event"""
+  # in case no event is provided, find it using eventNumber
+  if event is None:
+    if eventNumber is None:
+      print "DumpEventInfo Error: either pass an event or an event number"
+      return
+    # find event based on run and event
+    if os.path.isdir(path):
+      dirList=os.listdir(path)
+      files=[]
+      for fname in dirList:
+        files.append(path+fname)
+    elif os.path.isfile(path):
+      files=[path]
+    else:
+      files=[]
+    events = AnalysisEvent(files)
+    EventSelection.prepareAnalysisEvent(events)
+    DumpEventInfo(events[eventNumber])
+    return
+  # run the producers when we want to print the outcome, and mute unneeded collections
+  for product in eventDumpConfig.productsToPrint: getattr(event,product)
+  for collection in eventDumpConfig.collectionsToHide: event.removeCollection(collection)
+  # Now, we can go on with the printing.
+  print event
+
+if __name__=="__main__":
+  import sys
+  DumpEventInfo(eventNumber=int(sys.argv[1]), path=sys.argv[2])
+
Index: python/DelphesAnalysis/DumpEventList.py
===================================================================
--- python/DelphesAnalysis/DumpEventList.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/DumpEventList.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,35 @@
+#!/usr/bin/env python 
+import os
+import ROOT
+import Delphes
+from AnalysisEvent import AnalysisEvent
+import EventSelection
+
+def DumpEventList(category, path="./", output="eventlist.txt"):
+  """Dump a list of events in a given category"""
+  # input
+  if os.path.isdir(path):
+    dirList=os.listdir(path)
+    files=[]
+    for fname in dirList:
+      files.append(path+fname)
+  elif os.path.isfile(path):
+    files=[path]
+  else:
+    files=[]
+  # events
+  events = AnalysisEvent(files)
+  # output
+  event_list = open(output,"w")
+  # collections and producers used in the analysis
+  EventSelection.prepareAnalysisEvent(events)
+  for event in events:
+    # check category
+    if EventSelection.isInCategory(category, event.category):
+      # print
+      print >> event_list , "Event", event.event()
+
+if __name__=="__main__":
+  import sys
+  DumpEventList(int(sys.argv[1]), path=sys.argv[2], output="eventlist.txt")
+
Index: python/DelphesAnalysis/EventSelection.py
===================================================================
--- python/DelphesAnalysis/EventSelection.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/EventSelection.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,61 @@
+from CPconfig import configuration
+from importlib import import_module
+
+# A proper implementation of the event selection will at least 
+# define the following:
+# - categoryNames (list of names of categories)
+# - eventCategory (category data producer)
+# - isInCategory (category definition, from category data)
+
+# the list of category names
+categoryNames = [ ]
+
+def eventCategory(event):
+  """Check analysis requirements for various steps
+     and return a tuple of data used to decide 
+     to what category an event belong """
+  raise NotImplementedError
+  return []
+
+def isInCategory(category, categoryData):
+  """Check if the event enters category X, given the tuple computed by eventCategory."""
+  raise NotImplementedError
+  return True
+
+# specific implementation of the "virtual methods" above
+EventSelectionImplementation = __import__(configuration.eventSelection)
+categoryNames  = EventSelectionImplementation.categoryNames
+eventCategory  = EventSelectionImplementation.eventCategory
+isInCategory   = EventSelectionImplementation.isInCategory
+
+# Functions below should not be touched in any implementation.
+
+def eventCategories():
+  """Return the number of catefories """
+  return len(categoryNames)
+
+def categoryName(category):
+  """Return the name of category x"""
+  if category<eventCategories() and category>=0: return categoryNames[category]
+  else: return "None"
+
+def categoriesHierarchy():
+  """Return a hierarchically organized list of categories"""
+  dct = {}
+  for i,item in enumerate(categoryNames):
+    p = dct
+    splitted = item.split('/')
+    for x in splitted[:-1]:
+      p = p.setdefault(x,{})
+    p[splitted[-1]] = i
+  return dct
+
+def prepareAnalysisEvent(event):
+  """Define collections and producers"""
+  for coll in configuration.eventCollections:
+    event.addCollection(coll.label,coll.collection)
+  for prod in configuration.eventProducers:
+    event.addProducer(prod.label,getattr(import_module(prod.module),prod.function),**prod.kwargs)
+  for weight in configuration.eventWeights:
+    events.addWeight(weight.label,getattr(import_module(weight.module),weight.classname)(**weight.kwargs))
+
Index: python/DelphesAnalysis/EventSelectionControlPlots.py
===================================================================
--- python/DelphesAnalysis/EventSelectionControlPlots.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/EventSelectionControlPlots.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,42 @@
+import ROOT
+import sys
+import os
+from AnalysisEvent import AnalysisEvent
+from BaseControlPlots import BaseControlPlots
+from CPconfig import configuration
+import EventSelection
+
+# Requirements:
+# fully implemented EventSelection
+# categories stored as event.category
+
+class EventSelectionControlPlots(BaseControlPlots):
+    """A class to create control plots for event selection"""
+
+    def __init__(self, dir=None, dataset=None, mode="plots"):
+      # create output file if needed. If no file is given, it means it is delegated
+      BaseControlPlots.__init__(self, dir=dir, purpose="eventSelection", dataset=dataset, mode=mode)
+      self.eventCategories = EventSelection.eventCategories()
+ 
+    def beginJob(self):
+      # declare histograms
+      self.add("event","Event number",1000,0,100000)
+      self.add("category","event category",self.eventCategories+1,0,self.eventCategories+1)
+
+    def process(self, event):
+      """EventSelectionControlPlots"""
+      result = { }
+      ## event category
+      categoryData = event.category
+      result["category"] = [ ]
+      for category in range(self.eventCategories):
+        if EventSelection.isInCategory(category, categoryData):
+          result["category"].append(category)
+      result["event"] = event.event()
+      return result
+
+if __name__=="__main__":
+  import sys
+  from BaseControlPlots import runTest
+  runTest(sys.argv[1], EventSelectionControlPlots())
+
Index: python/DelphesAnalysis/__init__.py
===================================================================
--- python/DelphesAnalysis/__init__.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
+++ python/DelphesAnalysis/__init__.py	(revision ae1d49f529a4c42441725d12c11a5ef5774c264b)
@@ -0,0 +1,6 @@
+__all__ = [ "Delphes", 
+            "AnalysisEvent", 
+            "BaseControlPlots", "BaseWeightClass", 
+            "EventSelection", 
+            "EventSelectionControlPlots", 
+            "ControlPlots", "DumpEventInfo" ]
