feat: support iguana::Algorithm::Run calls

Clas12Banks/clas12reader.cpp

@@ -56,6 +56,7 @@ namespace clas12 {
       addBank(additionalBank);
 
     _applyQA=other._applyQA;
+    _readEventUserAction=other._readEventUserAction;
   }
 
   void clas12reader::initReader(){
@@ -461,8 +462,9 @@ namespace clas12 {
       _event.getStructure(*ibank.get());
     }
 
-
-    return true;
+    // now that we have read all the banks, call user's custom read action;
+    // return its return value, since it's the last thing `readEvent()` does
+    return _readEventUserAction(this);
   }
   ////////////////////////////////////////////////////////
   ///initialise next event from the reader
@@ -644,18 +646,22 @@ namespace clas12 {
   }
 
    ////////////////////////////////////////////////////////
-  ///Filter and return detParticles by given PID
-  std::vector<region_part_ptr> clas12reader::getByID(int id){
-    return container_filter(_detParticles, [id](region_part_ptr dr)
-			    {return dr->getPid()==id;});
+  ///Filter and return detParticles
+  std::vector<region_part_ptr> clas12reader::getByID(int id, bool applyBankFilter){
+    return container_filter(_detParticles, [id,applyBankFilter](region_part_ptr dr)
+			    {return dr->getPid()==id && (!applyBankFilter || dr->isAllowed());});
+  }
+  std::vector<region_part_ptr> clas12reader::getByRegion(int ir, bool applyBankFilter){
+    return container_filter(_detParticles, [ir,applyBankFilter](region_part_ptr dr)
+			    {return dr->getRegion()==ir && (!applyBankFilter || dr->isAllowed());});
   }
-  std::vector<region_part_ptr> clas12reader::getByRegion(int ir){
-    return container_filter(_detParticles, [ir](region_part_ptr dr)
-			    {return dr->getRegion()==ir;});
+  std::vector<region_part_ptr> clas12reader::getByCharge(int ch, bool applyBankFilter){
+    return container_filter(_detParticles, [ch,applyBankFilter](region_part_ptr dr)
+			    {return dr->par()->getCharge()==ch && (!applyBankFilter || dr->isAllowed());});
   }
-  std::vector<region_part_ptr> clas12reader::getByCharge(int ch){
-    return container_filter(_detParticles, [ch](region_part_ptr dr)
-			    {return dr->par()->getCharge()==ch;});
+  std::vector<region_part_ptr> clas12reader::getDetParticles(bool applyBankFilter){
+    return container_filter(_detParticles, [applyBankFilter](region_part_ptr dr)
+			    {return applyBankFilter ? dr->isAllowed() : true;});
   }
 
   ////////////////////////////////////////////////////////////////

Clas12Banks/clas12reader.h

@@ -44,6 +44,7 @@
 #include <vector>
 #include <string>
 #include <iostream>
+#include <functional>
 
 #ifdef RCDB_MYSQL
    #include "rcdb_reader.h"
@@ -114,7 +115,31 @@ namespace clas12 {
       _rbands.push_back(std::move(reg));
      }
 
-
+    // bank accessors
+    // - the method name convention is `getBankName` where `BankName` is the bank name without the colons (`::`)
+    // - the `class` tag is needed disambiguate between classes and `clas12reader` methods with the same name
+    class particle&      getRECParticle()      const { return *_bparts;     }
+    class ftbparticle&   getRECFTParticle()    const { return *_bftbparts;  }
+    class helonline&     getHELonline()        const { return *_bhelonline; }
+    class helflip&       getHELflip()          const { return *_bhelflip;   }
+    class runconfig&     getRUNconfig()        const { return *_brunconfig; }
+    class event&         getRECEvent()         const { return *_bevent;     }
+    class ftbevent&      getRECFTEvent()       const { return *_bftbevent;  }
+    class vtp&           getRAWvtp()           const { return *_bvtp;       }
+    class vertdoca&      getRECVertDoca()      const { return *_bvertdoca;  }
+    class calorimeter&   getRECCalorimeter()   const { return *_bcal;       }
+    class scintillator&  getRECScintillator()  const { return *_bscint;     }
+    class tracker&       getRECTrack()         const { return *_btrck;      }
+    class covmatrix&     getRECCovMat()        const { return *_bcovmat;    }
+    class utracker&      getRECUTrack()        const { return *_butrck;     }
+    class traj&          getRECTraj()          const { return *_btraj;      }
+    class cherenkov&     getRECCherenkov()     const { return *_bcher;      }
+    class rich&          getRICHParticle()     const { return *_brich;      }
+    class forwardtagger& getRECForwardTagger() const { return *_bft;        }
+    class mcparticle&    getMCLund()           const { return *_bmcparts;   }
+    class mcevent&       getMCEvent()          const { return *_bmcevent;   }
+
+    // bank pointer accessors
     helonline_ptr helonline() const{return _bhelonline.get();};
     helflip_ptr helflip() const{return _bhelflip.get();};
     runconfig_ptr runconfig() const{return _brunconfig.get();};
@@ -170,12 +195,39 @@ namespace clas12 {
     }
 
     /////////////////////////////////
-
+    // region particle accessors
+
+    /// @returns the number of region particles; note that this does not take into account any `hipo::bank` filters (_e.g._, from iguana)
+    int getNParticles() const  noexcept{return _detParticles.size();}
+
+    /// @returns a reference to the full list of region particles
+    /// @see getDetParticles(bool) if you need a `hipo::bank` filter applied
     std::vector<region_part_ptr>& getDetParticles(){return _detParticles;}
+
+    /// @returns a pointer to the full list of region particles
+    /// @see getDetParticles(bool) if you need a `hipo::bank` filter applied
     std::vector<region_part_ptr>* getDetParticlesPtr(){return &_detParticles;}
-    std::vector<region_part_ptr> getByID(int id);
-    std::vector<region_part_ptr> getByRegion(int ir);
-    std::vector<region_part_ptr> getByCharge(int ch);
+
+    /// @returns a _copy_ of the full list of region particles, with or without applying a `hipo::bank` filter
+    /// @param applyBankFilter if `true`, apply any `hipo::bank` filters (_e.g._, from iguana)
+    std::vector<region_part_ptr> getDetParticles(bool applyBankFilter=false);
+
+    /// @returns a subset of region particles, filtered by PDG
+    /// @param id the PDG
+    /// @param applyBankFilter if `true`, apply any `hipo::bank` filters (_e.g._, from iguana)
+    std::vector<region_part_ptr> getByID(int id, bool applyBankFilter=false);
+
+    /// @returns a subset of region particles, filtered by region
+    /// @param ir the region
+    /// @param applyBankFilter if `true`, apply any `hipo::bank` filters (_e.g._, from iguana)
+    std::vector<region_part_ptr> getByRegion(int ir, bool applyBankFilter=false);
+
+    /// @returns a subset of region particles, filtered by charge
+    /// @param ch the charge
+    /// @param applyBankFilter if `true`, apply any `hipo::bank` filters (_e.g._, from iguana)
+    std::vector<region_part_ptr> getByCharge(int ch, bool applyBankFilter=false);
+
+    /////////////////////////////////
 
     const std::vector<short>& preCheckPids();
     const std::vector<short>& preCheckPidsOrCharge();
@@ -195,7 +247,6 @@ namespace clas12 {
 
     void useFTBased(){_useFTBased=true;}
 
-    int getNParticles() const  noexcept{return _detParticles.size();}
     const std::vector<short> &getPids() const  noexcept{return _pids;}
 
     bool checkTriggerBit(uint k){
@@ -257,6 +308,15 @@ namespace clas12 {
       _verbose=level;
       _reader.setVerbose(level);
     }
+
+    /// Set a "read action", a custom lambda function that is executed for every event within `readEvent()`,
+    /// which is called by methods like `clas12reader::next()` and `HipoChain::Next()`
+    /// @param readEventUserAction lambda function, where its argument is a pointer to an instance of this `clas12reader` class,
+    /// and its `bool` return value controls whether to proceed with the event or not
+    void SetReadAction(std::function<bool(clas12reader*)> readEventUserAction) {
+      _readEventUserAction = readEventUserAction;
+    }
+
 
   protected:
 
@@ -345,7 +405,10 @@ namespace clas12 {
     bool _useFTBased{false};
     bool _isOpen{false};
     std::vector<std::string> _addBankNames;
-
+
+    // user-definable lambda, called in `readEvent()`; the default is a no-op returning `true`
+    std::function<bool(clas12reader*)> _readEventUserAction = [](clas12reader* r) {return true;};
+
     ///////////////////////////////DB
   private:
     int _runNo{0};

Clas12Banks/region_particle.h

@@ -61,6 +61,9 @@ namespace clas12 {
     /// i.e. how the detector banks relate to that region
     virtual bool sort(){
       _pentry=_parts->getEntry();
+      //check if allowed by particle-bank filters
+      _allowed=bankAllowsRow(_pentry, _parts);
+      _allowed_ftb=bankAllowsRow(_pentry, _ftbparts);
       //check for covariance matrix
       if(_covmat)_pcmat=_covmat->getIndex(_pentry);
       if(_mcpart)_pmc=_mcpart->match_to(_pentry);
@@ -158,9 +161,27 @@ namespace clas12 {
     double getMCPhiDiff() const {return getPhi()-mc()->getPhi();}
     double getMCPDiff() const {return getP()-mc()->getP();}
 
+    /// @returns true if this particle is "allowed", _e.g._, by a HIPO bank filter,
+    /// which can be applied by an iguana algorithm
+    bool const& isAllowed() const {
+      if(_ftbparts==nullptr) return _allowed;
+      return _useFTBPid*_ftbparts->getRows() ? _allowed_ftb : _allowed;
+    }
+
     //if(_parts->getCharge())
   protected:
 
+    /// @returns `true` if `bank` contains row `row` in its `hipo::bank::rowlist`
+    /// @param row the row to look up
+    /// @param bank the bank to check
+    bool const bankAllowsRow(int const& row, hipo::bank const* bank) const {
+      // NOTE: `hipo::bank::getRowList()` returns the list of rows which are _allowed_ by `hipo::bank::filter()`; for example,
+      // `hipo::bank::filter` is used by iguana fiducial cut algorithms to select particles from `REC::Particle`
+      if(bank!=nullptr && bank.getRows()>0)
+        return std::find(bank->getRowList().begin(), bank->getRowList().end(), row) != bank->getRowList().end();
+      return false;
+    }
+
     par_ptr _parts={nullptr};
     ftbpar_ptr _ftbparts={nullptr};
     covmat_ptr _covmat={nullptr};
@@ -183,6 +204,10 @@ namespace clas12 {
     short _pcmat=-1;
     short _region=-1;
     short _useFTBPid=0;
+
+    // filter status (from `hipo::bank::filter`)
+    bool _allowed{true}; // allowed by `_parts`, i.e., `REC::Particle`
+    bool _allowed_ftb{true}; // allowed by `_ftbparts`, i.e., `RECFT::Particle`
   };
   //pointer "typedef"
   using region_part_ptr=clas12::region_particle*;

README.md

@@ -578,63 +578,10 @@ Note the event number is just its position in the file, not the DST RUN::Config:
 ## Ex 11 Iguana interface
 
 To run iguana routines in clas12root you should first set the environment to point to an
-installed version of iguana, by setting the IGUANA variable to the INSTALLATION directory.
-
-There are 2 methods of using Iguana. For speed and simplicity both use just the iguana
-action function and will not operate on the underlying banks structures.
-In the first method, which can use any iguan algorithm, you just directly create the iguana
-algorithm and run it yourself. See Ex11_Iguana_example_01_bank_rows.C
-
-The second method is though a higher level interface which simplifies the usage, but is not
-garuanteed to be able to use all algorithms. In these cases you may revert to the first method
-to apply any additional algorithms. The interface for this is in the $CLAS12ROOT/iguana/ directory
-and an example is available at Ex11_Iguana_MomentumCorrection.C.
-
-These example use a HipoChain and you will need to set file paths yourself etc.
-Using clas12root:region_particles means that information for
-alogithms is readily available and so we can just pass these particle objects
-into the action functions.
-
-Highlighting the iguana parts for method 2 :
-
-      //clas12root-iguana interface
-      clas12root::Iguana ig;
-
-      //choose some algorithms to apply
-      ig.GetTransformers().Use("clas12::MomentumCorrection");
-      ig.GetFilters().Use("clas12::ZVertexFilter");
-      ig.GetCreators().Use("physics::InclusiveKinematics");
-
-      ig.SetOptionAll("log", "debug");
-      ig.Start();
-
-
-      ...
-      // get particles by type
-      // note we applied a filter to ensure size of all ==1
-      auto electron=c12->getByID(11)[0];
-      auto pip=c12->getByID(211)[0];
-      auto pim=c12->getByID(-211)[0];
-
-      ///Now do some iguana analysis!!
-
-      //filter on z-vertices of the particles
-      //note I can pass a vector of clas12root particles
-      if( !(ig.GetFilters().doZVertexFilter({electron,pip,pim})) ) {
-        continue;
-      }
-
-
-      //correct momentum and get 4-vectors
-      //I pass a vector of clas12root particles and LorentzVectors
-      ig.GetTransformers().doMomentumCorrections({electron,pip,pim},{&p4el,&p4pip,&p4pim});
-
-      //calculate inclusive kinematics
-      //use original bank info
-      auto kine = ig.GetCreators().doInclusiveKinematics(electron);
-      //use momentum corrected momentum
-      auto corrkine = ig.GetCreators().doInclusiveKinematics(p4el);
-
-      ...
-
+installed version of iguana, by setting the `IGUANA` variable to the INSTALLATION directory.
+Clas12root does _not_ depend on Iguana, but running `clas12root` will load its libraries for you if
+you have set the `IGUANA` environment variable.
 
+For usage of Iguana with `clas12root`, see the examples, such as
+- [`Ex11_Iguana.C`](/RunRoot/Ex11_Iguana.C)
+- [`Ex11_Iguana_example_01_bank_rows.C`](/RunRoot/Ex11_Iguana_example_01_bank_rows.C)