Changeset b7b836a in git for external/fastjet/contribs

Timestamp:

Jun 6, 2018, 10:05:10 PM (6 years ago)

Author:

Pavel Demin <pavel-demin@…>

Branches:

ImprovedOutputFile, Timing, dual_readout, llp, master

Children:

17d0ab8

Parents:

95e6b7a

Message:

update FastJet library to 3.3.1 and FastJet Contrib library to 1.036

Location:

external/fastjet/contribs

Files:

• Nsubjettiness/AUTHORS (modified) (1 diff)
• Nsubjettiness/AxesDefinition.hh (modified) (14 diffs)
• Nsubjettiness/ChangeLog (modified) (1 diff)
• Nsubjettiness/MeasureDefinition.cc (modified) (1 diff)
• Nsubjettiness/NEWS (modified) (1 diff)
• Nsubjettiness/VERSION (modified) (1 diff)
• RecursiveTools/AUTHORS (modified) (2 diffs)
• RecursiveTools/BottomUpSoftDrop.cc (added)
• RecursiveTools/BottomUpSoftDrop.hh (added)
• RecursiveTools/COPYING (modified) (1 diff)
• RecursiveTools/ChangeLog (modified) (1 diff)
• RecursiveTools/IteratedSoftDrop.cc (added)
• RecursiveTools/IteratedSoftDrop.hh (added)
• RecursiveTools/ModifiedMassDropTagger.hh (modified) (3 diffs)
• RecursiveTools/NEWS (modified) (1 diff)
• RecursiveTools/README (modified) (5 diffs)
• RecursiveTools/Recluster.hh (modified) (2 diffs)
• RecursiveTools/RecursiveSoftDrop.cc (added)
• RecursiveTools/RecursiveSoftDrop.hh (added)
• RecursiveTools/RecursiveSymmetryCutBase.cc (modified) (10 diffs)
• RecursiveTools/RecursiveSymmetryCutBase.hh (modified) (11 diffs)
• RecursiveTools/SoftDrop.cc (modified) (2 diffs)
• RecursiveTools/SoftDrop.hh (modified) (2 diffs)
• RecursiveTools/TODO (modified) (2 diffs)
• RecursiveTools/VERSION (modified) (1 diff)
• ValenciaPlugin/Makefile (deleted)

external/fastjet/contribs/Nsubjettiness/AUTHORS

@@ -39 +39 @@
    Iain W. Stewart, Frank J. Tackmann, Jesse Thaler,
    Christopher K. Vermilion, and Thomas F. Wilkason.
-   arXiv:1508.01516.
+   JHEP 1511:072 (2015), arXiv:1508.01516.
 
    Resolving Boosted Jets with XCone.
    Jesse Thaler and Thomas F. Wilkason.
-   arXiv:1508.01518.
+   JHEP 1512:051 (2015), arXiv:1508.01518.
 
 -------------------------------------------------------------------------------

external/fastjet/contribs/Nsubjettiness/AxesDefinition.hh

@@ -5 +5 @@
 //  Jesse Thaler, Ken Van Tilburg, Christopher K. Vermilion, and TJ Wilkason
 //
-//  $Id: AxesDefinition.hh 833 2015-07-23 14:35:23Z jthaler $
+//  $Id: AxesDefinition.hh 1130 2018-06-06 12:09:46Z jthaler $
 //----------------------------------------------------------------------
 // This file is part of FastJet contrib.
@@ -508 +508 @@
 /// \brief Wrapper for jet definitions (for memory management)
 ///
-/// This class was introduced to avoid issue of a FastJet bug when using genKT clustering
-/// Now using this for all AxesDefinition with a manual recombiner to use the delete_recombiner_when_unused function
+/// This class is used by all AxesDefinition with a manual recombiner to
+/// ensure that the delete_recombiner_when_unused function is always called
 ///------------------------------------------------------------------------
 class JetDefinitionWrapper {
@@ -549 +549 @@
    : ExclusiveJetAxes(JetDefinitionWrapper(fastjet::kt_algorithm,
                                           fastjet::JetDefinition::max_allowable_R, //maximum jet radius constant
-                                          _recomb = new WinnerTakeAllRecombiner(), // Needs to be explicitly declared (this will be deleted by JetDefinitionWrapper)
+                                          new WinnerTakeAllRecombiner(), // Needs to be explicitly declared (this will be deleted by JetDefinitionWrapper)
                                           fastjet::Best).getJetDef()
                       ) {
@@ -571 +571 @@
    virtual WTA_KT_Axes* create() const {return new WTA_KT_Axes(*this);}
 
-private:
-   const WinnerTakeAllRecombiner *_recomb;  ///< Internal recombiner
-
 };
    
@@ -588 +585 @@
    : ExclusiveJetAxes(JetDefinitionWrapper(fastjet::cambridge_algorithm,
                                              fastjet::JetDefinition::max_allowable_R, //maximum jet radius constant
-                                             _recomb = new WinnerTakeAllRecombiner(), // Needs to be explicitly declared (this will be deleted by JetDefinitionWrapper)
+                                             new WinnerTakeAllRecombiner(), // Needs to be explicitly declared (this will be deleted by JetDefinitionWrapper)
                                              fastjet::Best).getJetDef()) {
     setNPass(NO_REFINING);
@@ -608 +605 @@
    /// For copying purposes
    virtual WTA_CA_Axes* create() const {return new WTA_CA_Axes(*this);}
-   
-private:
-   const WinnerTakeAllRecombiner *_recomb;  ///< Internal recombiner
-
+  
 };
 
@@ -674 +668 @@
                                             R0,
                                             p,
-                                            _recomb = new WinnerTakeAllRecombiner()
+                                            new WinnerTakeAllRecombiner()
                                             ).getJetDef()), _p(p), _R0(R0) {
       if (p < 0) throw Error("WTA_GenKT_Axes:  Currently only p >=0 is supported.");
@@ -702 +696 @@
    double _p;   ///< genkT power
    double _R0;  ///< jet radius
-   const WinnerTakeAllRecombiner *_recomb; ///< Internal recombiner
 };
    
@@ -717 +710 @@
    /// Constructor
    GenET_GenKT_Axes(double delta, double p, double R0 = fastjet::JetDefinition::max_allowable_R)
-   : ExclusiveJetAxes((JetDefinitionWrapper(fastjet::genkt_algorithm, R0, p, _recomb = new GeneralEtSchemeRecombiner(delta))).getJetDef() ),
+   : ExclusiveJetAxes((JetDefinitionWrapper(fastjet::genkt_algorithm, R0, p, new GeneralEtSchemeRecombiner(delta))).getJetDef() ),
     _delta(delta), _p(p), _R0(R0) {
        if (p < 0) throw Error("GenET_GenKT_Axes:  Currently only p >=0 is supported.");
@@ -750 +743 @@
    double _p;     ///< GenkT power
    double _R0;    ///< jet radius
-   const GeneralEtSchemeRecombiner *_recomb;   ///< Internal recombiner
 };
 
@@ -1201 +1193 @@
    /// Constructor
    Comb_WTA_GenKT_Axes(int nExtra, double p, double R0 = fastjet::JetDefinition::max_allowable_R)
-   : ExclusiveCombinatorialJetAxes((JetDefinitionWrapper(fastjet::genkt_algorithm, R0, p, _recomb = new WinnerTakeAllRecombiner())).getJetDef(), nExtra),
+   : ExclusiveCombinatorialJetAxes((JetDefinitionWrapper(fastjet::genkt_algorithm, R0, p, new WinnerTakeAllRecombiner())).getJetDef(), nExtra),
     _p(p), _R0(R0) {
        if (p < 0) throw Error("Comb_WTA_GenKT_Axes:  Currently only p >=0 is supported.");
@@ -1227 +1219 @@
    double _p;        ///< GenkT power
    double _R0;       ///< jet radius
-   const WinnerTakeAllRecombiner *_recomb;   ///< Internal recombiner
 };
    
@@ -1241 +1232 @@
    /// Constructor
    Comb_GenET_GenKT_Axes(int nExtra, double delta, double p, double R0 = fastjet::JetDefinition::max_allowable_R)
-   : ExclusiveCombinatorialJetAxes((JetDefinitionWrapper(fastjet::genkt_algorithm, R0, p, _recomb = new GeneralEtSchemeRecombiner(delta))).getJetDef(), nExtra),
+   : ExclusiveCombinatorialJetAxes((JetDefinitionWrapper(fastjet::genkt_algorithm, R0, p, new GeneralEtSchemeRecombiner(delta))).getJetDef(), nExtra),
     _delta(delta), _p(p), _R0(R0) {
        if (p < 0) throw Error("Comb_GenET_GenKT_Axes:  Currently only p >=0 is supported.");
@@ -1271 +1262 @@
    double _p;        ///< GenkT power
    double _R0;       ///< jet radius
-   const GeneralEtSchemeRecombiner *_recomb;  ///<  Internal recombiner
 };
    

external/fastjet/contribs/Nsubjettiness/ChangeLog

@@ +1 @@
+2018-07-06 <jthaler>
+   Updated comments in AxesDefinition.hh about role of JetDefinitionWrapper
+   Updated AUTHORS with JHEP publication information for XCone
+   Prepared VERSION and NEWS for 2.2.5 release
+2018-07-05 <jthaler>
+   Fixed bug in AxesDefinition.hh where _recomb was used before it was declared.
 2016-06-08 <jthaler>
    Fixed bug in MeasureDefinition.cc where axes were not completely defined,

external/fastjet/contribs/Nsubjettiness/MeasureDefinition.cc

@@ -5 +5 @@
 //  Jesse Thaler, Ken Van Tilburg, Christopher K. Vermilion, and TJ Wilkason
 //
-//  $Id: MeasureDefinition.cc 946 2016-06-14 19:11:27Z jthaler $
+//  $Id: MeasureDefinition.cc 943 2016-06-09 00:55:29Z jthaler $
 //----------------------------------------------------------------------
 // This file is part of FastJet contrib.

external/fastjet/contribs/Nsubjettiness/NEWS

@@ -32 +32 @@
     N-jettiness as a jet finder using the new ConicalGeometric measure.
 
+-- 2.2.5:  (June 6, 2018) Fixed bug involved undefined pointer for in
+           AxesDefinition (thanks Attila Krasznahorkay)
 -- 2.2.4:  (Jun 14, 2016) Fixed bug where multi-pass minimization could yield
            pathological axes (thanks Gregory Soyez) 

external/fastjet/contribs/Nsubjettiness/VERSION

@@ -1 @@
-2.2.4
+2.2.5

external/fastjet/contribs/RecursiveTools/AUTHORS

@@ -4 +4 @@
    Gregory Soyez <soyez@fastjet.fr>
    Jesse Thaler <jthaler@jthaler.net>
+   Kevin Zhou <knzhou@mit.edu>
+   Frederic Dreyer <fdreyer@mit.edu>
 
 The physics is based on:
@@ -16 +18 @@
    Andrew J. Larkoski, Simone Marzani, Gregory Soyez, and Jesse Thaler.
    JHEP 1405:146 (2014), arXiv:1402.2657
+
+   [IteratedSoftDrop]
+   Casimir Meets Poisson: Improved Quark/Gluon Discrimination with Counting Observables.
+   Christopher Frye, Andrew J. Larkoski, Jesse Thaler, Kevin Zhou.
+   JHEP 1709:083 (2017), arXiv:1704.06266
+
+   [RecursiveSoftDrop]
+   [BottomUpSoftDrop]
+   Recursive Soft Drop.
+   Frederic A. Dreyer, Lina Necib, Gregory Soyez, and Jesse Thaler
+   arXiv:1804.03657
+   

external/fastjet/contribs/RecursiveTools/COPYING

@@ -14 +14 @@
   arXiv:1307.0007
   arXiv:1402.2657
+  arXiv:1704.06266
 
 ======================================================================

external/fastjet/contribs/RecursiveTools/ChangeLog

@@ +1 @@
+2018-05-29  Jesse Thaler  <jthaler@jthaler.net>
+
+        * VERSION
+        * NEWS
+        Changed to 2.0.0-beta2, noted in news
+
+2018-04-21  Jesse Thaler  <jthaler@jthaler.net>
+
+        * AUTHORS: updated arxiv number for RecursiveSoftDrop
+
+2018-04-21  Gavin Salam  <gavin.salam@cern.ch>
+
+        * README: updated arxiv number for RecursiveSoftDrop & Bottom-up
+        Soft Drop.
+
+2018-04-04  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSoftDrop.cc (contrib):
+        fixed syntax of calls to structure_of<...>
+        (thanks to Attila Krasznahorkay)
+
+2018-01-24  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSoftDrop.cc:
+        for the (dafault) dynamical R0 implementation, the dynamical R0 is
+        evolved independently in each branch.
+
+2017-10-10  Jesse Thaler  <jthaler@jthaler.net>
+        * AUTHORS
+        Added mention of BUSD
+
+        * README
+        Some tweaks to the wording
+
+2017-10-11  Gregory Soyez  <soyez@fastjet.fr>
+
+        * IteratedSoftDrop.hh:
+        IteratedSoftDropInfo::size() and multiplicity() now return an
+        unsigned int instead of a double
+
+2017-10-10  Jesse Thaler  <jthaler@jthaler.net>
+        * AUTHORS
+        Updated journal reference for ISD
+
+        * example_isd.{cc,ref}:
+        Included soft drop multiplicity in example
+
+        * README
+        Added warning at top that documentation has not been updated
+
+        * VERSION
+        Changed to 2.0.0-beta1
+
+2017-10-10  Gregory Soyez  <soyez@fastjet.fr>
+
+        * example_isd.ref:
+        updated to reflect the bugfix below
+
+        * IteratedSoftDrop.cc:
+        fixed issue in description (was taking sqrt of -ve number when
+        there were no angular cut)
+
+        * NEWS:
+        drafted for RecursiveTools-2.0.0-beta1
+
+        * TODO:
+        updated list in view of a beta release
+
+        * example_isd.cc:
+        pointed to the right header for IteratedSoftDrop
+
+        * RecursiveSoftDrop.{hh,cc}:
+        * IteratedSoftDrop.{hh,cc}:
+        moved IteratedSoftDrop to its own file (replacing the old
+        implementation)
+
+        * example_advanced_usage.ref:
+        updated reference file following the fix below.
+
+2017-09-28  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSymmetryCutBase.cc:
+        when no substructure is found, keep the _symmetry, _delta_R and
+        _mu structure variables at -1.  This for example allows one to
+        trigger on substructure by checking if delta_R>=0.
+
+        Note: we could set it to 0 (as it was before) and trigger using
+        _delta_R>0 but there might be some genuine substructure exactly
+        collinear.
+
+2017-09-19  Gregory Soyez  <soyez@fastjet.fr>
+
+        * example_isd.ref:
+        updated to the latest version of the example
+
+2017-09-18  Gregory Soyez  <soyez@fastjet.fr>
+
+        * Makefile:
+        updating make check to use all the examples (failing on ISD as
+        expected, see below)
+
+        * example_bottomup_softdrop.cc:
+        fixed typo
+
+        * example_bottomup_softdrop.ref:
+        * example_recursive_softdrop.ref:
+        added reference output for this example
+
+        * RecursiveSoftDrop.{hh,cc}:
+        * RecursiveSymmetryCutBase.{cc,hh}:
+        moved the "all_prongs" method from the base structure t oa
+        standalone function in RecursiveSoftDrop.hh
+
+        In practice, this is irrelevant for mMDT and SD (since pieces()
+        gets the job done, and the substructure class does not have (as
+        is) reliable info to get the full structure)
+
+        * RecursiveSymmetryCutBase.cc:
+        revamped a series of requests for substructure info to better
+        handle possible recursion into deeper jet substructure.
+
+        * RecursiveSoftDrop.{hh,cc}:
+        updated "description" to reuse the info from the base class
+
+        * example_isd.cc:
+        updated to use the newer implementation of ISD. Checked that it
+        gives the same results as the former implementation.
+
+        Note: make check still needs fixing because the example now
+              computes a different set of angularities
+
+        * RecursiveSoftDrop.hh:
+        added a few helpers to IteratedSoftDropInfo ([] operator and size,
+        meaning it can be used as a vector<pair<double,double> >)
+
+        * RecursiveSymmetryCutBase.cc:
+        . fixed bugs in the calculation of the geometric distances for ee
+          coordinates
+        . fixed bug in the computation of the (zg,thetag) pairs [it was
+          returning the groomed ones instead of the ones passing the
+          condition]
+
+        * example_recursive_softdrop.cc:
+        set the R0 parameter to the original jet radius
+
+2017-09-13  Gregory Soyez  <soyez@fastjet.fr>
+
+        * example_recursive_softdrop.cc:
+        tied up a few comments and the code output
+
+        * RecursiveSymmetryCutBase.{hh,cc}:
+        removed the unneeded _is_composite
+
+        * RecursiveSoftDrop.cc:
+        fixed issue with "verbose" dropped info on branches with no
+        further substructure
+
+2017-09-11  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSoftDrop.{hh,cc}:
+        have IteratedSoftDDrop returning by default an object of type
+        IteratedSoftDropInfo; added several helpers
+
+2017-09-08  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSoftDrop.{hh,cc}:
+        updated IteratedSoftDrop to give it the flexibility of
+        RecursiveSoftDrop
+
+        * RecursiveSymmetryCutBase.hh:
+        fixed typo in comment
+
+        * example_mmdt_ee.cc:     *** ADDED ***
+        added an example to illustrat usage in ee collisions
+
+        * example_isd.cc:
+        * BottomUpSoftDrop.cc:
+        * IteratedSoftDrop.cc:
+        * RecursiveSoftDrop.cc:
+        Fixed compilation issues with FJ 3.0 (mostly the usage of features
+        introduced only in FJ3.1)
+
+        * RecursiveSymmetryCutBase.{hh,cc}:
+        used the internal Recluster class for FJ<3.1.0 and the FJ antive
+        one for FJ>=3.1.0
+
+        * BottomUpSoftDrop.{hh,cc}:
+        moved the implementation of global_grooming to the source file and
+        fixed a few trivial compilation errors
+
+2017-09-08  Frédéric Dreyer  <frederic.dreyer@gmail.com>
+
+        * BottomUpSoftDrop.hh:
+        added the global_grooming method to process full event
+
+2017-09-07  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSoftDrop.cc:
+        cleaned (mostly by removing older commented-out code)
+
+        * RecursiveSoftDrop.{hh,cc}:
+        * RecursiveSymmetryCutBase.{hh,cc}:
+        * SoftDrop.cc:
+
+        added support for ee coordinates.  For that, the symmetry measure
+        has to be set to either theta_E (which uses the 3-vector angle
+        theta) or to cos_theta_E which uses sqrt(2*[1-cos(theta)])
+
+        Accordingly, the recursion_choice can be set to larger_E to
+        recurse into the branch with the largest energy. The larger_m
+        mode, recorsing into the larger-mass branch is also possible but
+        not advised (for the same reason as the pp version).
+
+        * RecursiveSymmetryCutBase.{hh,cc}:
+        switched to the Recluster class provided with FastJet. ASlso
+        included the recluster description to RecursiveSymmetryCutBase
+        when it is user-defined.
+
+2017-09-06  Gregory Soyez  <soyez@fastjet.fr>
+
+        * BottomUpSoftDrop.{hh,cc}:
+        . renamed SoftDropStructure  -> BottomUpSoftDropStructure
+                  SoftDropRecombiner -> BottomUpSoftDropRecombiner
+                  SoftDropPlugin     -> BottomUpSoftDropPlugin
+        . moved 'description' to source file (instead of header)
+        . kept the "area" information when available (jets will just
+        appear as having a 0 area)
+        . added most class description (main class still missing)
+
+        * RecursiveSoftDrop.cc:
+        . put internal classes to handle CS history in an internal namespace
+        . replaced the "switch" in the mail loop by a series of if (allows
+          us a few simplificatins/cleaning)
+        . more uniform treatment of issues in the presence of an angular cut
+          (as part of the above reorganisation)
+
+        * example_advanced_usage.ref:
+        updated reference output file following the bugfix (missing
+        "grooming mode" initialisation in one of the SoftDrop ctors) on
+        2017-08-01
+
+        * RecursiveSymmetryCutBase.cc:
+        removed redundent code
+
+2017-08-10  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSoftDrop.cc:
+        fixed trivial typo in variable name
+
+>>>>>>> .r1071
+2017-08-04  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSoftDrop.cc:
+        do not impose an angular cut in IterativeSD if it is -ve.
+
+2017-08-01  Gregory Soyez  <soyez@fastjet.fr>
+
+        * example_recursive_softdrop.cc:
+        added a series of optional flags
+
+        * RecursiveSoftDrop.cc:
+        fixed a few issues with the fixed depth version
+
+        * RecursiveSymmetryCutBase.hh:
+        a jet is now considered as havig substructure if deltaR>0
+        (coherent with released version)
+
+        * SoftDrop.hh:
+        bugfix: set the "grooming mode" by default in all ctors
+        EDIT: caused issue with make check, fixed on 2017-09-069 (see
+              above)
+
+        * RecursiveSoftDrop.{hh,cc}:
+        added support for the "same depth" variant
+
+        * RecursiveSymmetryCutBase.cc:
+        also associate a RecursiveSymmetryCutBase::StructureType structure
+        to the result jet in case it is just a single particle (in
+        grooming mode)
+
+2017-07-31  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSymmetryCutBase.{hh,cc}:
+        added the option to pass an extra parameter to the symmetry cut
+        function
+
+        * RecursiveSymmetryCutBase.{hh,cc}:
+        * ModifiedMassDropTagger.hh
+        * SoftDrop.hh:
+        minor adaptions due to the above change + added a few methods to
+        query the class information (symmetry cut, beta, ...)
+
+        * RecursiveSoftDrop.{hh,cc}:
+        Added support for
+        - a dynamical R0
+        - recursing only in the hardest branch
+        - imposing a min deltaR cut
+
+        Added a tentative IterativeSoftDrop class
+
+2017-07-28  Gregory Soyez  <soyez@fastjet.fr>
+
+        * RecursiveSoftDrop.cc:
+        adapted to the latest changes in RecursiveSymmetryCutBase
+
+        * RecursiveSymmetryCutBase.{hh,cc}:
+        reorganised the output of the recursion step (recurse_one_step)
+        using an enum to make iot more readable (and to fix issues where
+        the dropped prong is actually 0, e.g. after subtraction)
+
+2017-07-26  Gregory Soyez  <soyez@fastjet.fr>
+
+        * example_recursive_softdrop.cc:   *** ADDED ***
+        added a draft example for the use of RecursiveSoftDrop
+
+        * RecursiveSoftDrop.{hh,cc}:       *** ADDED ***
+        added a first pass at an implementation of RecursiveSoftDrop.
+
+        This is based on Frederic's initial version but keeps the
+        branching structure of the jet. Some of the features, like
+        dynamical R0, direct access to the substructure or the same depth
+        variant, are still unsupported.
+
+        * SoftDrop.hh:
+        declared _beta, _symmetry_cut and _R0sqr as protected (was
+        private) so they ca n be used by RecursiveSoftDrop
+
+        * RecursiveSymmetryCutBase.{hh,cc}:
+        extracted from result() the part that performs one step of the
+        resursion (implemented as recurse_one_step()). This is repeatedly
+        called by result(). It has specific conventions to indicate
+        whether or not some substructure has been found or if one ran into
+        some issue.
+
+2017-04-25  Kevin Zhou  <knzhou@mit.edu>
+
+        * IteratedSoftDrop.hh
+          . Added Doxygen documentation
+
+        * RecursiveSymmetryCutBase.hh
+          . Added references to ISD
+
+2017-04-25  Jesse Thaler  <jthaler@jthaler.net>
+        * AUTHORS, COPYING:
+          . Added ISD arXiv number
+
+        * example_isd.{cc,ref}
+          . Added ISD arXiv number
+          . Changing z_cut to be optimal value (with Lambda = 1 GeV)
+          . Tweaked formatting
+
+        * IteratedSoftDrop.{hh,cc}
+          . Added ISD arXiv number
+          . Assert added if recluster does not return one jet.
+
+        * README
+          . Added ISD arXiv number and tweaked wording
+
+
+2017-04-20  Kevin Zhou  <knzhou@mit.edu>
+
+        * example_isd.{cc,ref}      ** ADDED **
+        * IteratedSoftDrop.{cc,hh}  ** ADDED **
+
+        * Makefile
+          . Added IteratedSoftDrop (ISD) as appropriate
+
+        * AUTHORS
+          . Added myself to author list
+          . Added placeholder for ISD paper
+
+        * COPYING
+          . Added placeholder for ISD paper
+
+        * README
+          . Added description of ISD
+
+        * TODO
+          . Added tasks to integrate ISD with other classes, e.g. SD, 
+            Recluster, and the future RecursiveSoftDrop (RSD)
+
+        * NEWS
+          . Added dummy for release of 1.1.0
+
+        * VERSION:
+          . Switched version number to 1.1.0-dev
+
+        * example_advanced_usage.cc:
+          . Added $Id$ tag, didn't change anything else
+
+2014-07-30  Gregory Soyez  <soyez@fastjet.fr>
+
+        * Recluster.hh: fixed the name of the #define for the header
+
 2014-07-09  Gregory Soyez  <soyez@fastjet.fr> + Jesse
 

external/fastjet/contribs/RecursiveTools/ModifiedMassDropTagger.hh

@@ -1 @@
-// $Id: ModifiedMassDropTagger.hh 688 2014-06-17 14:29:56Z jthaler $
+// $Id: ModifiedMassDropTagger.hh 1032 2017-07-31 14:20:03Z gsoyez $
 //
 // Copyright (c) 2014-, Gavin P. Salam
@@ -105 +105 @@
   virtual ~ModifiedMassDropTagger(){}
 
+  //----------------------------------------------------------------------
+  // access to class info
+  double symmetry_cut() const { return _symmetry_cut; }
+  
 protected:
 
@@ -110 +114 @@
   /// has no dependence on the subjet kinematics
   virtual double symmetry_cut_fn(const PseudoJet & /* p1 */, 
-                                 const PseudoJet & /* p2 */
+                                 const PseudoJet & /* p2 */,
+                                 void *extra_parameters = 0
                                  ) const {return _symmetry_cut;}
   virtual std::string symmetry_cut_description() const;

external/fastjet/contribs/RecursiveTools/NEWS

@@ +1 @@
+2018/05/31: release of version 2.0.0-beta2 with corrected syntax
+
+2017/10/10: release of version 2.0.0-beta1 including implementations of
+
+* RecursiveSoftDrop (see example_rsd.hh for usage)
+* IteratedSoftDrop (see example_isd.hh for usage)
+* e+e- version of the recursive tools (see example_mmdt_ee.hh for usage)
+* BottomUpSoftDrop (see example_bottomup_softdrop.cc for usage)
+
 2014/07/09: release of version 1.0.0 of RecursiveTools including
             ModifiedMassDropTagger and SoftDrop (as well as Recluster)

external/fastjet/contribs/RecursiveTools/README

@@ -17 +17 @@
   Marzani, Gregory Soyez, Jesse Thaler
 
+- RecursiveSoftDrop
+- BottomUpSoftDrop
+  This corresponds to arXiv:1804.03657 by Frederic Dreyer, Lina
+  Necib, Gregory Soyez and Jesse Thaler
+
+- IteratedSoftDrop
+  This corresponds to arXiv:1704.06266 by Christopher Frye, Andrew J.
+  Larkoski, Jesse Thaler, Kevin Zhou
+
 - Recluster
   A generic tool to recluster a given jet into subjets
-  Note: this is largely based on the Filter code in FastJet v3.0 and
-  ultimately, this tool will probably be moved into FastJet
+  Note: a Recluster class is available natively in FastJet since v3.1.
+        Users are therefore encouraged to use the FastJet version
+        rather than this one which is mostly provided for
+        compatibility of this contrib with older versions of FastJet.
 
 The interface for these tools is described in more detail below, with
@@ -74 +85 @@
 
 The SoftDrop procedure is very similar to mMDT, albeit with a 
-generalized symmetry condition:
+generalised symmetry condition:
  
    z > z_cut * (R / R0)^beta
@@ -108 +119 @@
 
 ------------------------------------------------------------------------
+RecursiveSoftDrop
+------------------------------------------------------------------------
+
+The RecursiveSoftDrop procedure applies the Soft Drop procedure N times 
+in a jet in order to find up to N+1 prongs.  N=0 makes no modification
+to the jet, and N=1 is equivalent to the original SoftDrop.
+
+Once one has more than one prong, one has to decide which will be
+declustered next.  At each step of the declustering procedure, one
+undoes the clustering which has the largest declustering angle
+(amongst all the branches that are searched for substructure). [see
+"set_fixed_depth" below for an alternative]
+
+Compared to SoftDrop, RecursiveSoftDrop takes an extra argument N
+specifying the number of times the SoftDrop procedure is recursively
+applied. Negative N means that the procedure is applied until no
+further substructure is found (i.e. corresponds to taking N=infinity).
+
+   double z_cut = 0.10;
+   double beta  = 2.0;
+   double R0    = 1.0; // this is the default value
+   int N        = -1;  
+   RecursiveSoftDrop rsd(z_cut, beta, N, R0);
+
+One then acts on a jet as
+
+   PseudoJet groomed_jet = rsd(jet)
+
+and get additional information via
+
+   groomed_jet.structure_of<RecursiveSoftDrop>()
+
+------------------------------------------------------------------------
+IteratedSoftDrop
+------------------------------------------------------------------------
+
+Iterated Soft Drop (ISD) is a repeated variant of SoftDrop.  After
+performing the Soft Drop procedure once, it logs the groomed symmetry
+factor, then recursively performs Soft Drop again on the harder
+branch.  This procedure is repeated down to an (optional) angular cut
+theta_cut, yielding a set of symmetry factors from which observables
+can be built.
+
+An IteratedSoftDrop tool can be created as follows:
+
+  double beta = -1.0;
+  double z_cut = 0.005;
+  double theta_cut = 0.0;
+  double R0 = 0.5;        // characteristic radius of jet algorithm
+  IteratedSoftDrop isd(beta, z_cut, double theta_cut, R0);
+
+By default, ISD applied on a jet gives a result of type
+IteratedSoftDropInfo that can then be probed to obtain physical
+observables
+
+  IteratedSoftDropInfo isd_info = isd(jet);
+
+  unsigned int multiplicity = isd_info.multiplicity();
+  double kappa = 1.0;     // changes angular scale of ISD angularity
+  double isd_width = isd_info.angularity(kappa);
+  vector<pair<double,double> > zg_thetags = isd_info.all_zg_thetag();
+  vector<pair<double,double> > zg_thetags = isd_info();
+  for (unsigned int i=0; i< isd_info.size(); ++i){
+    cout << "(zg, theta_g)_" << i << " = "
+         << isd_info[i].first << " " <<  isd_info[i].second << endl;
+  }
+
+Alternatively, one can directly get the multiplicity, angularity, and
+(zg,thetag) pairs from the IteratedSoftDrop class, at the expense of
+re-running the declustering procedure:
+
+  unsigned int multiplicity = isd.multiplicity(jet);
+  double isd_width = isd.angularity(jet, 1.0);
+  vector<pair<double,double> > zg_thetags = isd.all_zg_thetag(jet);
+
+
+Note: the iterative declustering procedure is the same as what one
+  would obtain with RecursiveSoftDrop with an (optional) angular cut
+  and recursing only in the hardest branch [see the "Changing
+  behaviour" section below for details], except that it returns some
+  information about the jet instead of a modified jet as RSD does.
+
+
+------------------------------------------------------------------------
+BottomUpSoftDrop
+------------------------------------------------------------------------
+
+This is a bottom-up version of the RecursiveSoftDrop procedure, in a
+similar way as Pruning can be seen as a bottom-up version of Trimming.
+
+In practice, the jet is reclustered and at each step of the clustering
+one checks the SoftDrop condition
+
+   z > z_cut * (R / R0)^beta
+
+If the condition is met, the pair is recombined. If the condition is
+not met, only the hardest of the two objects is kept for further
+clustering and the softest is rejected.
+
+------------------------------------------------------------------------
 Recluster
 ------------------------------------------------------------------------
+
+  *** NOTE: this is provided only for backwards compatibility ***
+  *** with FastJet <3.1. For FastJet >=3.1, the native        ***
+  *** fastjet::Recluster is used instead                      ***
 
 The Recluster class allows the constituents of a jet to be reclustered
 with a different recursive clustering algorithm.  This is used
-internally in the mMDT/SoftDrop code in order to recluster the jet using
-the CA algorithm.  This is achieved via
+internally in the mMDT/SoftDrop/RecursiveSoftDrop/IteratedSoftDrop
+code in order to recluster the jet using the CA algorithm.  This is
+achieved via
 
   Recluster ca_reclusterer(cambridge_algorithm,
@@ -123 +239 @@
 delete_self_when_unused.
 
-
 ------------------------------------------------------------------------
 Changing behaviour
 ------------------------------------------------------------------------
 
-The behaviour of the ModifiedMassDropTagger and SoftDrop classes can
-be tweaked using the following options:
-
-SymmetryMeasure = {scalar_z,vector_z,y}  [constructor argument]
+The behaviour of the all the tools provided here
+(ModifiedMassDropTagger, SoftDrop, RecursiveSoftDrop and
+IteratedSoftDrop) can be tweaked using the following options:
+
+SymmetryMeasure = {scalar_z, vector_z, y, theta_E, cos_theta_E}
+  [constructor argument]
   : The definition of the energy sharing between subjets, with 0
-    corresponding to the most asymmetric
-
-RecursionChoice = {larger_pt,larger_mt,larger_m}  [constructor argument]
+    corresponding to the most asymmetric.
+    . scalar_z = min(pt1,pt2)/(pt1+pt2)   [default]
+    . vector_z = min(pt1,pt2)/pt_{1+2}
+    . y = min(pt1^2,pt2^2)/m_{12}^2  (original y from MDT)
+    . theta_E     = min(E1,E2)/(E1+E2),
+           with angular measure theta_{12}^2
+    . cos_theta_E = min(E1,E2)/(E1+E2),
+           with angular measure 2[1-cos(theta_{12})]
+    The last two variants are meant for use in e+e- collisions,
+    together with the "larger_E" recursion choice (see below)
+
+RecursionChoice = {larger_pt, larger_mt, larger_m, larger_E}
+  [constructor argument]
   : The path to recurse through the tree after the symmetry condition
-    fails
+    fails. Options refer to transverse momentum (pt), transverse mass
+    (mt=sqrt(pt^2+m^2), mass (m) or energy (E). the latter is meant
+    for use in e+e- collisions
 
 mu_cut   [constructor argument]
   : An optional mass drop condition
 
-set_subtractor(subtractor*) [or subtracter as a constructor argument]
+set_subtractor(subtractor*) [or subtractor as a constructor argument]
   : provide a subtractor. When a subtractor is supplied, the
     kinematic constraints are applied on subtracted 4-vectors. In
@@ -165 +294 @@
     and SoftDrop defaults to grooming mode.
 
+set_verbose_structure(bool)
+  : when set to true, additional information will be stored in the jet
+    structure. This includes in particular values of symmetry,
+    delta_R, and mu of dropped branches
+
+For the specific case of RecursiveSoftDrop, additional tweaking is
+possible via the following methods
+
+set_fixed_depth_mode(bool)
+  : when this is true, RSD will recurse (N times) into all the
+    branches found during the previous iteration [instead of recursing
+    through the largest declustering angle until N prongs have been
+    found]. This yields at most 2^N prong. For infinite N, the two
+    options are equivalent.
+
+set_dynamical_R0(bool)
+  : By default the angles in the SD condition are normalised to the
+    parameter R0. With "dynamical R0", RSD will dynamically adjust R0
+    to be the angle between the two prongs found during the previous
+    iteration.
+
+set_hardest_branch_only(bool)
+  : When substructure is found, only recurse into the hardest of the
+    two branches for further substructure search. This uses the class
+    RecursionChoice.
+
+set_min_deltaR_squared(double):
+  : set a minimal angle (squared) at which we stop the declustering
+    procedure.  This cut is ineffective for negative values of the
+    argument.
+
 ------------------------------------------------------------------------
 Technical Details

external/fastjet/contribs/RecursiveTools/Recluster.hh

@@ -1 @@
-#ifndef __FASTJET_TOOLS_RECLUSTER_HH__
-#define __FASTJET_TOOLS_RECLUSTER_HH__
+#ifndef __FASTJET_CONTRIB_TOOLS_RECLUSTER_HH__
+#define __FASTJET_CONTRIB_TOOLS_RECLUSTER_HH__
 
-// $Id: Recluster.hh 700 2014-07-07 12:50:05Z gsoyez $
+// $Id: Recluster.hh 723 2014-07-30 09:11:01Z gsoyez $
 //
 // Copyright (c) 2014-, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
@@ -181 +181 @@
 FASTJET_END_NAMESPACE      // defined in fastjet/internal/base.hh
 
-#endif   // __FASTJET_TOOLS_RECLUSTER_HH__
+#endif   // __FASTJET_CONTRIB_TOOLS_RECLUSTER_HH__

external/fastjet/contribs/RecursiveTools/RecursiveSymmetryCutBase.cc

@@ -1 @@
-// $Id: RecursiveSymmetryCutBase.cc 700 2014-07-07 12:50:05Z gsoyez $
+// $Id: RecursiveSymmetryCutBase.cc 1080 2017-09-28 07:51:37Z gsoyez $
 //
 // Copyright (c) 2014-, Gavin P. Salam, Gregory Soyez, Jesse Thaler
@@ -23 +23 @@
 #include "fastjet/JetDefinition.hh"
 #include "fastjet/ClusterSequenceAreaBase.hh"
+#include <cassert>
 #include <algorithm> 
 #include <cstdlib> 
@@ -42 +43 @@
 PseudoJet RecursiveSymmetryCutBase::result(const PseudoJet & jet) const {
   // construct the input jet (by default, recluster with C/A)
-
   if (! jet.has_constituents()){
     throw Error("RecursiveSymmetryCutBase can only be applied to jets with constituents");
   }
 
-  PseudoJet j = 
-    _do_reclustering 
-      ? _recluster ? (*_recluster)(jet)
-                   : Recluster(cambridge_algorithm, JetDefinition::max_allowable_R)(jet)
-      : jet;
-    
-  // issue a warning if the jet is not obtained through a C/A
-  // clustering
-  // if ((! j.has_associated_cluster_sequence()) ||
-  //     (j.validated_cs()->jet_def().jet_algorithm() != cambridge_algorithm))
-  //  _nonca_warning.warn("RecursiveSymmetryCutBase is designed to be applied on jets from a Cambridge/Aachen clustering; use it with other algorithms at your own risk.");
-
+  PseudoJet j = _recluster_if_needed(jet);
+
+  // sanity check: the jet must have a valid CS
   if (! j.has_valid_cluster_sequence()){
     throw Error("RecursiveSymmetryCutBase can only be applied to jets associated to a (valid) cluster sequence");
   }
 
+  // check that area information is there in case we have a subtractor
+  // GS: do we really need this since subtraction may not require areas?
   if (_subtractor) {
     const ClusterSequenceAreaBase * csab = 
@@ -75 +68 @@
     subjet = (*_subtractor)(subjet);
   }
-
-  bool use_mu_cut = (_mu_cut != numeric_limits<double>::infinity());
 
   // variables for tracking what will happen
@@ -86 +77 @@
   std::vector<double> dropped_symmetry;
   std::vector<double> dropped_mu;
+
+  double sym, mu2;
   
   // now recurse into the jet's structure
-  while (subjet.has_parents(piece1, piece2)) {
-    
-    // first sanity check: 
-    // - zero or negative pts are not allowed for the input subjet
-    // - zero or negative masses are not allowed for configurations
-    //   in which the mass will effectively appear in a denominator
-    //   (The masses will be checked later)
-    if (subjet.pt2() <= 0) return PseudoJet();
-
-    if (_subtractor) {
-      piece1 = (*_subtractor)(piece1);
-      piece2 = (*_subtractor)(piece2);
-    }
-    
-    // determine the symmetry parameter
-    double sym;
-
-    if        (_symmetry_measure == y) {
-      // the original d_{ij}/m^2 choice from MDT
-      // first make sure the mass is sensible
-      if (subjet.m2() <= 0) {
-        _negative_mass_warning.warn("RecursiveSymmetryCutBase: cannot calculate y, because (sub)jet mass is negative; bailing out");
-        return _result_no_substructure(subjet); //TBC: do we return the hardest parent? A NULL PseudoJet?
-      }
-      sym = piece1.kt_distance(piece2) / subjet.m2();
-
-    } else if (_symmetry_measure == vector_z) {
-      // min(pt1, pt2)/(pt), where the denominator is a vector sum
-      // of the two subjets
-      sym = min(piece1.pt(), piece2.pt()) / subjet.pt();
-
-    } else if (_symmetry_measure == scalar_z) {
-      // min(pt1, pt2)/(pt1+pt2), where the denominator is a scalar sum
-      // of the two subjets
-      double pt1 = piece1.pt();
-      double pt2 = piece2.pt();
-      // make sure denominator is non-zero
-      sym = pt1 + pt2;
-      if (sym == 0) return PseudoJet(); 
-      sym = min(pt1, pt2) / sym;
-
-    } else {
-      throw Error ("Unrecognized choice of symmetry_measure");
-    }
-
-    // determine the symmetry cut
-    // (This function is specified in the derived classes)
-    double this_symmetry_cut = symmetry_cut_fn(piece1, piece2);
-
-    // and make a first tagging decision based on symmetry cut
-    bool tagged = (sym > this_symmetry_cut);
-
-    // if tagged based on symmetry cut, then check the mu cut (if relevant)
-    // and update the tagging decision. Calculate mu^2 regardless, for cases
-    // of users not cutting on mu2, but still interested in its value.
-    double mu2 = max(piece1.m2(), piece2.m2())/subjet.m2();
-    if (tagged && use_mu_cut) {
-      // first a sanity check -- mu2 won't be sensible if the subjet mass 
-      // is negative, so we can't then trust the mu cut - bail out
-      if (subjet.m2() <= 0) {
-        _negative_mass_warning.warn("RecursiveSymmetryCutBase: cannot trust mu, because (sub)jet mass is negative; bailing out");
-        return PseudoJet();
-      }
-      if (mu2 > 1) _mu2_gt1_warning.warn("RecursiveSymmetryCutBase encountered mu^2 value > 1");
-      if (mu2 > pow(_mu_cut,2)) tagged = false;
-    }
-
-
-    // if we've tagged the splitting, return the jet with its substructure
-    if (tagged) {
-      // record relevant information
-      StructureType * structure = new StructureType(subjet);
-      structure->_symmetry = sym;
-      structure->_mu       = (mu2 >= 0) ? sqrt(mu2) : -sqrt(-mu2);
-      structure->_delta_R  = piece1.delta_R(piece2);
-      if (_verbose_structure) {
-        structure->_has_verbose = true;
-        structure->_dropped_symmetry = dropped_symmetry;
-        structure->_dropped_mu = dropped_mu;
-        structure->_dropped_delta_R = dropped_delta_R;
-      }
-      subjet.set_structure_shared_ptr(SharedPtr<PseudoJetStructureBase>(structure));
-      return subjet;
-    }
+  RecursionStatus status;
+  while ((status=recurse_one_step(subjet, piece1, piece2, sym, mu2)) != recursion_success) {
+    // start with sanity checks:
+    if ((status == recursion_issue) || (status == recursion_no_parents)) {
+      // we should return piece1 by our convention for recurse_one_step
+      PseudoJet result;
+      if (status == recursion_issue){
+        result = piece1;
+        if (_verbose) cout << "reached end; returning null jet " << endl;
+      } else {
+        result = _result_no_substructure(piece1);
+        if (_verbose) cout << "no parents found; returning last PJ or empty jet" << endl;
+      }
+      
+      if (result != 0) {
+        // if in grooming mode, add dummy structure information
+        StructureType * structure = new StructureType(result);
+        // structure->_symmetry = 0.0;
+        // structure->_mu       = 0.0;
+        // structure->_delta_R  = 0.0;
+        if (_verbose_structure) { // still want to store verbose information about dropped branches
+          structure->_has_verbose = true;
+          structure->_dropped_symmetry = dropped_symmetry;
+          structure->_dropped_mu = dropped_mu;
+          structure->_dropped_delta_R = dropped_delta_R;
+        }
+        result.set_structure_shared_ptr(SharedPtr<PseudoJetStructureBase>(structure));
+      }
+  
+      return result;
+    }
+
+    assert(status == recursion_dropped);
     
     // if desired, store information about dropped branches before recursing
@@ -179 +121 @@
       dropped_mu.push_back((mu2 >= 0) ? sqrt(mu2) : -sqrt(-mu2));
     }
-    
-    // otherwise continue unclustering, allowing for the different
-    // ways of choosing which parent to look into
-    int choice;
-    if        (_recursion_choice == larger_mt) {
-      choice = piece1.mt2() > piece2.mt2() ? 1 : 2;
-
-    } else if (_recursion_choice == larger_pt) {
-      choice = piece1.pt2() > piece2.pt2() ? 1 : 2;
-
-    } else if (_recursion_choice == larger_m)  {
-      choice = piece1.m2()  > piece2.m2()  ? 1 : 2;
-
-    } else {
-      throw Error ("Unrecognized value for recursion_choice");
-    }    
-    if (_verbose) cout << "choice is " << choice << endl;;
-    subjet = (choice == 1) ? piece1 : piece2;
-  } // (subjet.has_parents(...))
-
-  if (_verbose) cout << "reached end; returning null jet " << endl;
-  
-  // decide on tagging versus grooming mode here
-  PseudoJet result = _result_no_substructure(subjet);
-  
-  if (result != 0) {
-    // if in grooming mode, add dummy structure information
-    StructureType * structure = new StructureType(result);
-    structure->_symmetry = 0.0;
-    structure->_mu       = 0.0;
-    structure->_delta_R  = 0.0;
-    if (_verbose_structure) { // still want to store verbose information about dropped branches
-      structure->_has_verbose = true;
-      structure->_dropped_symmetry = dropped_symmetry;
-      structure->_dropped_mu = dropped_mu;
-      structure->_dropped_delta_R = dropped_delta_R;
-    }
-    result.set_structure_shared_ptr(SharedPtr<PseudoJetStructureBase>(structure));
-  }
-  
-  return result;
-}
-
-
+
+    subjet = piece1;
+  }
+  
+
+  // we've tagged the splitting, return the jet with its substructure
+  StructureType * structure = new StructureType(subjet);
+  structure->_symmetry = sym;
+  structure->_mu       = (mu2 >= 0) ? sqrt(mu2) : -sqrt(-mu2);
+  structure->_delta_R  = sqrt(squared_geometric_distance(piece1, piece2));
+  if (_verbose_structure) {
+    structure->_has_verbose = true;
+    structure->_dropped_symmetry = dropped_symmetry;
+    structure->_dropped_mu = dropped_mu;
+    structure->_dropped_delta_R = dropped_delta_R;
+  }
+  subjet.set_structure_shared_ptr(SharedPtr<PseudoJetStructureBase>(structure));
+  return subjet;
+}
+
+
+  
+//----------------------------------------------------------------------
+// the method below is the one actually performing one step of the
+// recursion.
+//
+// It returns a status code (defined above)
+//
+// In case of success, all the information is filled
+// In case of "no parents", piee1 is the same subjet
+// In case of trouble, piece2 will be a 0 PJ and piece1 is the PJ we
+//   should return (either 0 itself if the issue was critical, or
+//   non-wero in case of a minor issue just causing the recursion to
+//   stop)
+RecursiveSymmetryCutBase::RecursionStatus
+      RecursiveSymmetryCutBase::recurse_one_step(const PseudoJet & subjet,
+                                                 PseudoJet &piece1, PseudoJet &piece2,
+                                                 double &sym, double &mu2,
+                                                 void *extra_parameters) const {
+  if (!subjet.has_parents(piece1, piece2)){
+    piece1 = subjet;    
+    piece2 = PseudoJet();
+    return recursion_no_parents;
+  }
+
+  // first sanity check: 
+  // - zero or negative pts are not allowed for the input subjet
+  // - zero or negative masses are not allowed for configurations
+  //   in which the mass will effectively appear in a denominator
+  //   (The masses will be checked later)
+  if (subjet.pt2() <= 0){ // this is a critical problem, return an empty PJ
+    piece1 = piece2 = PseudoJet();
+    return recursion_issue;
+  }
+
+  if (_subtractor) {
+    piece1 = (*_subtractor)(piece1);
+    piece2 = (*_subtractor)(piece2);
+  }
+    
+  // determine the symmetry parameter
+  if        (_symmetry_measure == y) {
+    // the original d_{ij}/m^2 choice from MDT
+    // first make sure the mass is sensible
+    if (subjet.m2() <= 0) { 
+      _negative_mass_warning.warn("RecursiveSymmetryCutBase: cannot calculate y, because (sub)jet mass is negative; bailing out");
+      // since rounding errors can give -ve masses, be a it more
+      // tolerant and consider that no substructure has been found
+      piece1 = _result_no_substructure(subjet);
+      piece2 = PseudoJet();
+      return recursion_issue;
+    }
+    sym = piece1.kt_distance(piece2) / subjet.m2();
+    
+  } else if (_symmetry_measure == vector_z) {
+    // min(pt1, pt2)/(pt), where the denominator is a vector sum
+    // of the two subjets
+    sym = min(piece1.pt(), piece2.pt()) / subjet.pt();    
+  } else if (_symmetry_measure == scalar_z) {
+    // min(pt1, pt2)/(pt1+pt2), where the denominator is a scalar sum
+    // of the two subjets
+    double pt1 = piece1.pt();
+    double pt2 = piece2.pt();
+    // make sure denominator is non-zero
+    sym = pt1 + pt2;
+    if (sym == 0){ // this is a critical problem, return an empty PJ
+      piece1 = piece2 = PseudoJet();
+      return recursion_issue;
+    }
+    sym = min(pt1, pt2) / sym;
+  } else if ((_symmetry_measure == theta_E) || (_symmetry_measure == cos_theta_E)){
+    // min(E1, E2)/(E1+E2)
+    double E1 = piece1.E();
+    double E2 = piece2.E();
+    // make sure denominator is non-zero
+    sym = E1 + E2;
+    if (sym == 0){ // this is a critical problem, return an empty PJ
+      piece1 = piece2 = PseudoJet();
+      return recursion_issue;
+    }
+    sym = min(E1, E2) / sym;
+  } else {
+    throw Error ("Unrecognized choice of symmetry_measure");
+  }
+  
+  // determine the symmetry cut
+  // (This function is specified in the derived classes)
+  double this_symmetry_cut = symmetry_cut_fn(piece1, piece2, extra_parameters);
+  
+  // and make a first tagging decision based on symmetry cut
+  bool tagged = (sym > this_symmetry_cut);
+
+  // if tagged based on symmetry cut, then check the mu cut (if relevant)
+  // and update the tagging decision. Calculate mu^2 regardless, for cases
+  // of users not cutting on mu2, but still interested in its value.
+  bool use_mu_cut = (_mu_cut != numeric_limits<double>::infinity());
+  mu2 = max(piece1.m2(), piece2.m2())/subjet.m2();
+  if (tagged && use_mu_cut) {
+    // first a sanity check -- mu2 won't be sensible if the subjet mass 
+    // is negative, so we can't then trust the mu cut - bail out
+    if (subjet.m2() <= 0) {
+      _negative_mass_warning.warn("RecursiveSymmetryCutBase: cannot trust mu, because (sub)jet mass is negative; bailing out");
+      piece1 = piece2 = PseudoJet();
+      return recursion_issue;
+    }
+    if (mu2 > 1) _mu2_gt1_warning.warn("RecursiveSymmetryCutBase encountered mu^2 value > 1");
+    if (mu2 > pow(_mu_cut,2)) tagged = false;
+  }
+
+  // we'll continue unclustering, allowing for the different
+  // ways of choosing which parent to look into
+  if        (_recursion_choice == larger_pt) {
+    if (piece1.pt2() < piece2.pt2()) std::swap(piece1, piece2);    
+  } else if (_recursion_choice == larger_mt) {
+    if (piece1.mt2() < piece2.mt2()) std::swap(piece1, piece2);    
+  } else if (_recursion_choice == larger_m)  {
+    if (piece1.m2() < piece2.m2()) std::swap(piece1, piece2);    
+  } else if (_recursion_choice == larger_E)  {
+    if (piece1.E()  < piece2.E())  std::swap(piece1, piece2);    
+  } else {
+    throw Error ("Unrecognized value for recursion_choice");
+  }    
+
+  return tagged ? recursion_success : recursion_dropped;
+}
+
+  
 //----------------------------------------------------------------------
 string RecursiveSymmetryCutBase::description() const {
@@ -235 +280 @@
   case vector_z:
     ostr << "vector_z"; break;
+  case theta_E:
+    ostr << "theta_E"; break;
+  case cos_theta_E:
+    ostr << "cos_theta_E"; break;
   default:
     cerr << "failed to interpret symmetry_measure" << endl; exit(-1);
@@ -254 +303 @@
   case larger_m:
     ostr << "mass"; break;
+  case larger_E:
+    ostr << "energy"; break;
   default:
     cerr << "failed to interpret recursion_choice" << endl; exit(-1);
@@ -259 +310 @@
 
   if (_subtractor) {
-    ostr << " and subtractor: " << _subtractor->description();
+    ostr << ", subtractor: " << _subtractor->description();
     if (_input_jet_is_subtracted) {ostr << " (input jet is assumed already subtracted)";}
   }
+
+  if (_recluster) {
+    ostr << " and reclustering using " << _recluster->description();
+  }
+  
   return ostr.str();
 }
 
+//----------------------------------------------------------------------
+// helper for handling the reclustering
+PseudoJet RecursiveSymmetryCutBase::_recluster_if_needed(const PseudoJet &jet) const{
+  if (! _do_reclustering) return jet;
+  if (_recluster) return (*_recluster)(jet);
+  if (is_ee()){
+#if FASTJET_VERSION_NUMBER >= 30100
+    return Recluster(JetDefinition(ee_genkt_algorithm, JetDefinition::max_allowable_R, 0.0), true)(jet);
+#else
+    return Recluster(JetDefinition(ee_genkt_algorithm, JetDefinition::max_allowable_R, 0.0))(jet);
+#endif
+  }
+
+  return Recluster(cambridge_algorithm, JetDefinition::max_allowable_R)(jet);
+}  
+  
+//----------------------------------------------------------------------
 // decide what to return when no substructure has been found
+double RecursiveSymmetryCutBase::squared_geometric_distance(const PseudoJet &j1,
+                                                            const PseudoJet &j2) const{
+  if (_symmetry_measure == theta_E){
+    double dot_3d = j1.px()*j2.px() + j1.py()*j2.py() + j1.pz()*j2.pz();
+    double cos_theta = max(-1.0,min(1.0, dot_3d/sqrt(j1.modp2()*j2.modp2())));
+    double theta = acos(cos_theta);
+    return theta*theta;
+  } else if (_symmetry_measure == cos_theta_E){
+    double dot_3d = j1.px()*j2.px() + j1.py()*j2.py() + j1.pz()*j2.pz();
+    return max(0.0, 2*(1-dot_3d/sqrt(j1.modp2()*j2.modp2())));
+  }
+
+  return j1.squared_distance(j2);
+}
+
+//----------------------------------------------------------------------
 PseudoJet RecursiveSymmetryCutBase::_result_no_substructure(const PseudoJet &last_parent) const{
   if (_grooming_mode){
@@ -277 +366 @@
 
 
+//========================================================================
+// implementation of the details of the structure
+  
+// the number of dropped subjets
+int RecursiveSymmetryCutBase::StructureType::dropped_count(bool global) const {
+  check_verbose("dropped_count()");
+
+  // if this jet has no substructure, just return an empty vector
+  if (!has_substructure()) return _dropped_delta_R.size();
+
+  // deal with the non-global case
+  if (!global) return _dropped_delta_R.size();
+
+  // for the global case, we've unfolded the recursion (likely more
+  // efficient as it requires less copying)
+  unsigned int count = 0;
+  vector<const RecursiveSymmetryCutBase::StructureType*> to_parse;
+  to_parse.push_back(this);
+
+  unsigned int i_parse = 0;
+  while (i_parse<to_parse.size()){
+    const RecursiveSymmetryCutBase::StructureType *current = to_parse[i_parse];
+    count += current->_dropped_delta_R.size();
+
+    // check if we need to recurse deeper in the substructure
+    //
+    // we can have 2 situations here for the underlying structure (the
+    // one we've wrapped around):
+    //  - it's of the clustering type
+    //  - it's a composite jet
+    // only in the 2nd case do we have to recurse deeper
+    const CompositeJetStructure *css = dynamic_cast<const CompositeJetStructure*>(current->_structure.get());
+    if (css == 0){ ++i_parse; continue; }
+    
+    vector<PseudoJet> prongs = css->pieces(PseudoJet()); // argument irrelevant
+    assert(prongs.size() == 2);
+    for (unsigned int i_prong=0; i_prong<2; ++i_prong){
+      if (prongs[i_prong].has_structure_of<RecursiveSymmetryCutBase>()){
+        RecursiveSymmetryCutBase::StructureType* prong_structure
+          = (RecursiveSymmetryCutBase::StructureType*) prongs[i_prong].structure_ptr();
+        if (prong_structure->has_substructure())
+          to_parse.push_back(prong_structure);
+      }
+    }
+
+    ++i_parse;
+  }
+  return count;
+}
+
+// the delta_R of all the dropped subjets
+vector<double> RecursiveSymmetryCutBase::StructureType::dropped_delta_R(bool global) const {
+  check_verbose("dropped_delta_R()");
+
+  // if this jet has no substructure, just return an empty vector
+  if (!has_substructure()) return vector<double>();
+
+  // deal with the non-global case
+  if (!global) return _dropped_delta_R;
+
+  // for the global case, we've unfolded the recursion (likely more
+  // efficient as it requires less copying)
+  vector<double> all_dropped;
+  vector<const RecursiveSymmetryCutBase::StructureType*> to_parse;
+  to_parse.push_back(this);
+
+  unsigned int i_parse = 0;
+  while (i_parse<to_parse.size()){
+    const RecursiveSymmetryCutBase::StructureType *current = to_parse[i_parse];
+    all_dropped.insert(all_dropped.end(), current->_dropped_delta_R.begin(), current->_dropped_delta_R.end());
+
+    // check if we need to recurse deeper in the substructure
+    //
+    // we can have 2 situations here for the underlying structure (the
+    // one we've wrapped around):
+    //  - it's of the clustering type
+    //  - it's a composite jet
+    // only in the 2nd case do we have to recurse deeper
+    const CompositeJetStructure *css = dynamic_cast<const CompositeJetStructure*>(current->_structure.get());
+    if (css == 0){ ++i_parse; continue; }
+    
+    vector<PseudoJet> prongs = css->pieces(PseudoJet()); // argument irrelevant
+    assert(prongs.size() == 2);
+    for (unsigned int i_prong=0; i_prong<2; ++i_prong){
+      if (prongs[i_prong].has_structure_of<RecursiveSymmetryCutBase>()){
+        RecursiveSymmetryCutBase::StructureType* prong_structure
+          = (RecursiveSymmetryCutBase::StructureType*) prongs[i_prong].structure_ptr();
+        if (prong_structure->has_substructure())
+          to_parse.push_back(prong_structure);
+      }
+    }
+
+    ++i_parse;
+  }
+  return all_dropped;
+}
+
+// the symmetry of all the dropped subjets
+vector<double> RecursiveSymmetryCutBase::StructureType::dropped_symmetry(bool global) const {
+  check_verbose("dropped_symmetry()");
+
+  // if this jet has no substructure, just return an empty vector
+  if (!has_substructure()) return vector<double>();
+
+  // deal with the non-global case
+  if (!global) return _dropped_symmetry;
+
+  // for the global case, we've unfolded the recursion (likely more
+  // efficient as it requires less copying)
+  vector<double> all_dropped;
+  vector<const RecursiveSymmetryCutBase::StructureType*> to_parse;
+  to_parse.push_back(this);
+
+  unsigned int i_parse = 0;
+  while (i_parse<to_parse.size()){
+    const RecursiveSymmetryCutBase::StructureType *current = to_parse[i_parse];
+    all_dropped.insert(all_dropped.end(), current->_dropped_symmetry.begin(), current->_dropped_symmetry.end());
+
+    // check if we need to recurse deeper in the substructure
+    //
+    // we can have 2 situations here for the underlying structure (the
+    // one we've wrapped around):
+    //  - it's of the clustering type
+    //  - it's a composite jet
+    // only in the 2nd case do we have to recurse deeper
+    const CompositeJetStructure *css = dynamic_cast<const CompositeJetStructure*>(current->_structure.get());
+    if (css == 0){ ++i_parse; continue; }
+    
+    vector<PseudoJet> prongs = css->pieces(PseudoJet()); // argument irrelevant
+    assert(prongs.size() == 2);
+    for (unsigned int i_prong=0; i_prong<2; ++i_prong){
+      if (prongs[i_prong].has_structure_of<RecursiveSymmetryCutBase>()){
+        RecursiveSymmetryCutBase::StructureType* prong_structure
+          = (RecursiveSymmetryCutBase::StructureType*) prongs[i_prong].structure_ptr();
+        if (prong_structure->has_substructure())
+          to_parse.push_back(prong_structure);
+      }
+    }
+
+    ++i_parse;
+  }
+  return all_dropped;
+}
+
+// the mu of all the dropped subjets
+vector<double> RecursiveSymmetryCutBase::StructureType::dropped_mu(bool global) const {
+  check_verbose("dropped_mu()");
+
+  // if this jet has no substructure, just return an empty vector
+  if (!has_substructure()) return vector<double>();
+
+  // deal with the non-global case
+  if (!global) return _dropped_mu;
+
+  // for the global case, we've unfolded the recursion (likely more
+  // efficient as it requires less copying)
+  vector<double> all_dropped;
+  vector<const RecursiveSymmetryCutBase::StructureType*> to_parse;
+  to_parse.push_back(this);
+
+  unsigned int i_parse = 0;
+  while (i_parse<to_parse.size()){
+    const RecursiveSymmetryCutBase::StructureType *current = to_parse[i_parse];
+    all_dropped.insert(all_dropped.end(), current->_dropped_mu.begin(), current->_dropped_mu.end());
+
+    // check if we need to recurse deeper in the substructure
+    //
+    // we can have 2 situations here for the underlying structure (the
+    // one we've wrapped around):
+    //  - it's of the clustering type
+    //  - it's a composite jet
+    // only in the 2nd case do we have to recurse deeper
+    const CompositeJetStructure *css = dynamic_cast<const CompositeJetStructure*>(current->_structure.get());
+    if (css == 0){ ++i_parse; continue; }
+    
+    vector<PseudoJet> prongs = css->pieces(PseudoJet()); // argument irrelevant
+    assert(prongs.size() == 2);
+    for (unsigned int i_prong=0; i_prong<2; ++i_prong){
+      if (prongs[i_prong].has_structure_of<RecursiveSymmetryCutBase>()){
+        RecursiveSymmetryCutBase::StructureType* prong_structure
+          = (RecursiveSymmetryCutBase::StructureType*) prongs[i_prong].structure_ptr();
+        if (prong_structure->has_substructure())
+          to_parse.push_back(prong_structure);
+      }
+    }
+
+    ++i_parse;
+  }
+  return all_dropped;
+}
+
+// the maximum of the symmetry over the dropped subjets
+double RecursiveSymmetryCutBase::StructureType::max_dropped_symmetry(bool global) const {
+  check_verbose("max_dropped_symmetry()");
+
+  // if there is no substructure, just exit
+  if (!has_substructure()){ return 0.0; }
+
+  // local value of the max_dropped_symmetry
+  double local_max = (_dropped_symmetry.size() == 0)
+    ? 0.0 : *max_element(_dropped_symmetry.begin(),_dropped_symmetry.end());
+
+  // recurse down the structure if instructed to do so
+  if (global){
+    // we can have 2 situations here for the underlying structure (the
+    // one we've wrapped around):
+    //  - it's of the clustering type
+    //  - it's a composite jet
+    // only in the 2nd case do we have to recurse deeper
+    const CompositeJetStructure *css = dynamic_cast<const CompositeJetStructure*>(_structure.get());
+    if (css == 0) return local_max;
+    
+    vector<PseudoJet> prongs = css->pieces(PseudoJet()); // argument irrelevant
+    assert(prongs.size() == 2);
+    for (unsigned int i_prong=0; i_prong<2; ++i_prong){
+      // check if the prong has further substructure
+      if (prongs[i_prong].has_structure_of<RecursiveSymmetryCutBase>()){
+        RecursiveSymmetryCutBase::StructureType* prong_structure
+          = (RecursiveSymmetryCutBase::StructureType*) prongs[i_prong].structure_ptr();
+        local_max = max(local_max, prong_structure->max_dropped_symmetry(true));
+      }
+    }
+  }
+
+  return local_max;
+}
+
+//------------------------------------------------------------------------
+// helper class to sort by decreasing thetag
+class SortRecursiveSoftDropStructureZgThetagPair{
+public:
+  bool operator()(const pair<double, double> &p1, const pair<double, double> &p2) const{
+    return p1.second > p2.second;
+  }
+};
+//------------------------------------------------------------------------
+
+// the (zg,thetag) pairs of all the splitting that were found and passed the SD condition
+vector<pair<double,double> > RecursiveSymmetryCutBase::StructureType::sorted_zg_and_thetag() const {
+  //check_verbose("sorted_zg_and_thetag()");
+
+  // if this jet has no substructure, just return an empty vector
+  if (!has_substructure()) return vector<pair<double,double> >();
+  
+  // otherwise fill a vector with all the prongs (no specific ordering)
+  vector<pair<double,double> > all;
+  vector<const RecursiveSymmetryCutBase::StructureType*> to_parse;
+  to_parse.push_back(this);
+
+  unsigned int i_parse = 0;
+  while (i_parse<to_parse.size()){
+    const RecursiveSymmetryCutBase::StructureType *current = to_parse[i_parse];
+    all.push_back(pair<double,double>(current->_symmetry, current->_delta_R));
+
+    vector<PseudoJet> prongs = current->pieces(PseudoJet());
+    assert(prongs.size() == 2);
+    for (unsigned int i_prong=0; i_prong<2; ++i_prong){
+      if (prongs[i_prong].has_structure_of<RecursiveSymmetryCutBase>()){
+        RecursiveSymmetryCutBase::StructureType* prong_structure
+          = (RecursiveSymmetryCutBase::StructureType*) prongs[i_prong].structure_ptr();
+        if (prong_structure->has_substructure())
+          to_parse.push_back(prong_structure);
+      }
+    }
+
+    ++i_parse;
+  }
+
+  sort(all.begin(), all.end(), SortRecursiveSoftDropStructureZgThetagPair());
+  return all;
+}
+
 } // namespace contrib
 

external/fastjet/contribs/RecursiveTools/RecursiveSymmetryCutBase.hh

@@ -1 @@
-// $Id: RecursiveSymmetryCutBase.hh 700 2014-07-07 12:50:05Z gsoyez $
+// $Id: RecursiveSymmetryCutBase.hh 1074 2017-09-18 15:15:20Z gsoyez $
 //
 // Copyright (c) 2014-, Gavin P. Salam, Gregory Soyez, Jesse Thaler
@@ -24 +24 @@
 
 #include <limits>
+#include <cassert>
 #include <fastjet/internal/base.hh>
 #include "fastjet/tools/Transformer.hh"
 #include "fastjet/WrappedStructure.hh"
-
+#include "fastjet/CompositeJetStructure.hh"
+
+#include "fastjet/config.h"
+
+// we'll use the native FJ class for reculstering if available
+#if FASTJET_VERSION_NUMBER >= 30100
+#include "fastjet/tools/Recluster.hh"
+#else
 #include "Recluster.hh"
+#endif
 
 /** \mainpage RecursiveTools contrib 
 
-    The aims RecursiveTools contrib to provide a set of tools for
-    recursive investigation of the substructure of jets.
-
-    Currently it includes:
+    The RecursiveTools contrib provides a set of tools for
+    recursive investigation jet substructure. Currently it includes:
     - fastjet::contrib::ModifiedMassDropTagger
     - fastjet::contrib::SoftDrop
-    - the two above classes derive from fastjet::contrib::RecursiveSymmetryCutBase
-    - fastjet::contrib::Recluster provides a reclustering transformer
-    - example*.cc provides a usage examples
+    - fastjet::contrib::RecursiveSymmetryCutBase (from which the above two classes derive)
+    - fastjet::contrib::IteratedSoftDropSymmetryFactors (defines ISD procedure)
+    - fastjet::contrib::IteratedSoftDropMultiplicity (defines a useful observable using ISD)  
+    - example*.cc provides usage examples
  */
 
@@ -78 +86 @@
 
   /// an enum of the different (a)symmetry measures that can be used
-  enum SymmetryMeasure{scalar_z, ///< \f$ \min(p_{ti}, p_{tj})/(p_{ti} + p_{tj}) \f$
-                       vector_z, ///< \f$ \min(p_{ti}, p_{tj})/p_{t(i+j)} \f$
-                       y         ///  \f$ \min(p_{ti}^2,p_{tj}^2) \Delta R_{ij}^2 / m_{ij}^2 \f$
-                       
+  enum SymmetryMeasure{scalar_z,   ///< \f$ \min(p_{ti}, p_{tj})/(p_{ti} + p_{tj}) \f$
+                       vector_z,   ///< \f$ \min(p_{ti}, p_{tj})/p_{t(i+j)} \f$
+                       y,          ///< \f$ \min(p_{ti}^2,p_{tj}^2) \Delta R_{ij}^2 / m_{ij}^2 \f$
+                       theta_E,    ///< \f$ \min(E_i,E_j)/(E_i+E_j) \f$ with 3d angle (ee collisions)
+                       cos_theta_E ///< \f$ \min(E_i,E_j)/(E_i+E_j) \f$ with
+                                   ///  \f$ \sqrt{2[1-cos(theta)]}\f$ for angles (ee collisions)
   };
 
@@ -87 +97 @@
   enum RecursionChoice{larger_pt, ///< choose the subjet with larger \f$ p_t \f$
                        larger_mt, ///< choose the subjet with larger \f$ m_t \equiv (m^2+p_t^2)^{\frac12}] \f$
-                       larger_m   ///  choose the subjet with larger mass (deprecated)
+                       larger_m,  ///< choose the subjet with larger mass (deprecated)
+                       larger_E   ///< choose the subjet with larger energy (meant for ee collisions)
   }; 
 
@@ -119 +130 @@
   virtual ~RecursiveSymmetryCutBase(){}
 
+  // access to class info
+  //----------------------------------------------------------------------
+  SymmetryMeasure symmetry_measure() const { return _symmetry_measure; }
+  double mu_cut() const { return _mu_cut; }
+  RecursionChoice recursion_choice() const { return _recursion_choice; }
+
   // internal subtraction configuration
   //----------------------------------------------------------------------
@@ -147 +164 @@
   //----------------------------------------------------------------------
 
-  /// configure the reclustering prior to the SoftDrop de-clustering
+  /// configure the reclustering prior to the recursive de-clustering
   ///  \param do_reclustering   recluster the jet or not?
   ///  \param recluster         how to recluster the jet 
@@ -174 +191 @@
   /// in particular values of symmetry, delta_R, and mu of dropped branches
   void set_verbose_structure(bool enable=true) { _verbose_structure = enable; }
+  bool has_verbose_structure() const { return _verbose_structure; }
   
   
@@ -185 +203 @@
   /// subtracted jet.
   virtual PseudoJet result(const PseudoJet & j) const;
-
+  
   /// description of the tool
   virtual std::string description() const;
 
-  /// the type of the associated structure
-  //typedef RecursiveSymmetryCutBaseStructure StructureType;  
+  /// returns the gepometrical distance between the two particles
+  /// depending on the symmetry measure used
+  double squared_geometric_distance(const PseudoJet &j1,
+                                    const PseudoJet &j2) const;
+  
 
   class StructureType;
@@ -203 +224 @@
   /// apply for a given pair of subjets p1 and p2
   virtual double symmetry_cut_fn(const PseudoJet & /* p1 */, 
-                                 const PseudoJet & /* p2 */) const = 0;
+                                 const PseudoJet & /* p2 */,
+                                 void *extra_parameters = 0) const = 0;
   /// the associated dwescription
   virtual std::string symmetry_cut_description() const = 0;
+
+  //----------------------------------------------------------------------
+  /// this defines status codes when checking for substructure
+  enum RecursionStatus{
+    recursion_success=0,   //< found some substructure
+    recursion_dropped,     //< dropped softest prong; recursion continues
+    recursion_no_parents,  //< down to constituents; bottom of recursion
+    recursion_issue        //< something went wrong; recursion stops
+  };
+  
+  //----------------------------------------------------------------------
+  /// the method below is the one actually performing one step of the
+  /// recursion.
+  ///
+  /// It returns a status code (defined above)
+  ///
+  /// In case of success, all the information is filled
+  /// In case of "no parents", piee1 is the same subjet
+  /// In case of trouble, piece2 will be a 0 PJ and piece1 is the PJ we
+  ///   should return (either 0 itself if the issue was critical, or
+  ///   non-wero in case of a minor issue just causing the recursion to
+  ///   stop)
+  ///
+  /// The extra_parameter argument allows one to pass extra agruments
+  /// to the symmetry condition
+  RecursionStatus recurse_one_step(const PseudoJet & subjet,
+                                   PseudoJet &piece1, PseudoJet &piece2,
+                                   double &sym, double &mu2,
+                                   void *extra_parameters = 0) const;
+
+  //----------------------------------------------------------------------
+  /// helper for handling the reclustering
+  PseudoJet _recluster_if_needed(const PseudoJet &jet) const;
+  
+  //----------------------------------------------------------------------
+  // helpers for selecting between ee and pp cases
+  bool is_ee() const{
+    return ((_symmetry_measure==theta_E) || (_symmetry_measure==cos_theta_E));
+  }
 
 private:
@@ -238 +299 @@
 public:
   StructureType(const PseudoJet & j) :
-    WrappedStructure(j.structure_shared_ptr()),
+    WrappedStructure(j.structure_shared_ptr()), _delta_R(-1.0), _symmetry(-1.0), _mu(-1.0),
+    _is_composite(false), _has_verbose(false) // by default, do not store verbose structure
+  {}
+
+  /// construct a structure with
+  ///  - basic info inherited from the reference jet "j"
+  ///  - a given deltaR       for substructure
+  ///  - a given symmetry     for substructure
+  ///  - a given mu parameter for substructure
+  /// If j is a "copmposite jet", it means that it has further
+  /// substructure to potentially recurse  into
+  StructureType(const PseudoJet & j, double delta_R_in, double symmetry_in, double mu_in=-1.0) :
+    WrappedStructure(j.structure_shared_ptr()), _delta_R(delta_R_in), _symmetry(symmetry_in), _mu(mu_in),
+    _is_composite(dynamic_cast<const CompositeJetStructure*>(j.structure_ptr()) != NULL),
     _has_verbose(false) // by default, do not store verbose structure
   {}
   
   // information about kept branch
-  double delta_R()  const {return _delta_R;};
-  double symmetry() const {return _symmetry;};
-  double mu()       const {return _mu;};
+  double delta_R()  const {return _delta_R;}
+  double thetag()   const {return _delta_R;}  // alternative name
+  double symmetry() const {return _symmetry;}
+  double zg()       const {return _symmetry;} // alternative name
+  double mu()       const {return _mu;}
   
   // additional verbose information about dropped branches
   bool has_verbose() const { return _has_verbose;}
-  
-  // number of dropped branches
-  int dropped_count() const {
-    if (!_has_verbose) throw Error("RecursiveSymmetryCutBase::StructureType: Verbose structure must be turned on to get dropped_count() values.");
-    return _dropped_delta_R.size();
-  }
-  
-  // delta_R of dropped branches
-  std::vector<double> dropped_delta_R() const {
-    if (!_has_verbose) throw Error("RecursiveSymmetryCutBase::StructureType: Verbose structure must be turned on to get dropped_delta_R() values.");
-    return _dropped_delta_R;
-  }
-
-  // symmetry values of dropped branches
-  std::vector<double> dropped_symmetry() const {
-    if (!_has_verbose) throw Error("RecursiveSymmetryCutBase::StructureType: Verbose structure must be turned on to get dropped_symmetry() values.");
-    return _dropped_symmetry;
-  }
-
-  // mass drop values of dropped branches
-  std::vector<double> dropped_mu() const {
-    if (!_has_verbose) throw Error("RecursiveSymmetryCutBase::StructureType: Verbose structure must be turned on to get dropped_mu() values.");
-    return _dropped_mu;
-  }
-  
-  // maximum symmetry value dropped
-  double max_dropped_symmetry() const {
-    if (!_has_verbose) throw Error("RecursiveSymmetryCutBase::StructureType: Verbose structure must be turned on to get max_dropped_symmetry().");
-    if (_dropped_symmetry.size() == 0) return 0.0;
-    return *std::max_element(_dropped_symmetry.begin(),_dropped_symmetry.end());
-  }
-  
+  void set_verbose(bool value) { _has_verbose = value;}
+
+  /// returns true if the current jet has some substructure (i.e. has
+  /// been tagged by the resursion) or not
+  ///
+  /// Note that this should include deltaR==0 (e.g. a perfectly
+  /// collinear branching with SoftDrop)
+  bool has_substructure() const { return _delta_R>=0; }
+  
+  /// number of dropped branches
+  int dropped_count(bool global=true) const;
+  
+  /// delta_R of dropped branches
+  /// when "global" is set, recurse into possile further substructure
+  std::vector<double> dropped_delta_R(bool global=true) const;
+  void set_dropped_delta_R(const std::vector<double> &v) { _dropped_delta_R = v; }
+
+  /// symmetry values of dropped branches
+  std::vector<double> dropped_symmetry(bool global=true) const;
+  void set_dropped_symmetry(const std::vector<double> &v) { _dropped_symmetry = v; }
+
+  /// mass drop values of dropped branches
+  std::vector<double> dropped_mu(bool global=true) const;
+  void set_dropped_mu(const std::vector<double> &v) { _dropped_mu = v; }
+  
+  /// maximum symmetry value dropped
+  double max_dropped_symmetry(bool global=true) const;
+
+  /// (global) list of groomed away elements as zg and thetag
+  /// sorted from largest to smallest anlge
+  std::vector<std::pair<double,double> > sorted_zg_and_thetag() const;
+
 private:
   double _delta_R, _symmetry, _mu;
   friend class RecursiveSymmetryCutBase;
 
+  bool _is_composite;
+  
   // additional verbose information
   bool _has_verbose;
@@ -291 +369 @@
   std::vector<double> _dropped_symmetry;
   std::vector<double> _dropped_mu;
+
+  bool check_verbose(const std::string &what) const{
+    if (!_has_verbose){
+      throw Error("RecursiveSymmetryCutBase::StructureType: Verbose structure must be turned on to get "+what+".");
+      return false;
+    }
+    return true;
+  }
+    
   
 };

external/fastjet/contribs/RecursiveTools/SoftDrop.cc

@@ -1 @@
-// $Id: SoftDrop.cc 686 2014-06-14 03:25:09Z jthaler $
+// $Id: SoftDrop.cc 1059 2017-09-07 20:48:00Z gsoyez $
 //
 // Copyright (c) 2014-, Gregory Soyez, Jesse. Thaler
@@ -57 +57 @@
 //----------------------------------------------------------------------
 double SoftDrop::symmetry_cut_fn(const PseudoJet & p1, 
-                                 const PseudoJet & p2) const{
-  return _symmetry_cut * pow(p1.squared_distance(p2)/_R0sqr, 0.5*_beta);
+                                 const PseudoJet & p2,
+                                 void *optional_R0sqr_ptr) const{
+  double R0sqr = (optional_R0sqr_ptr == 0) ? _R0sqr : *((double*) optional_R0sqr_ptr);
+  return _symmetry_cut * pow(squared_geometric_distance(p1,p2)/R0sqr, 0.5*_beta);
 }
 

external/fastjet/contribs/RecursiveTools/SoftDrop.hh

@@ -1 @@
-// $Id: SoftDrop.hh 686 2014-06-14 03:25:09Z jthaler $
+// $Id: SoftDrop.hh 1034 2017-08-01 10:03:53Z gsoyez $
 //
 // Copyright (c) 2014-, Gregory Soyez, Jesse Thaler
@@ -112 +112 @@
      RecursiveSymmetryCutBase(symmetry_measure, mu_cut, recursion_choice, subtractor),
      _beta(beta), _symmetry_cut(symmetry_cut), _R0sqr(R0*R0)
-  {}
+  {
+    // change the default: use grooming mode
+    set_grooming_mode();
+  }
 
   /// default destructor
   virtual ~SoftDrop(){}
 
+  //----------------------------------------------------------------------
+  // access to class info
+  double beta()         const { return _beta; }
+  double symmetry_cut() const { return _symmetry_cut; }
+  double R0()           const { return sqrt(_R0sqr); }
+  
 protected:
 
   // Unlike MMDT, the SoftDrop symmetry_cut_fn depends on the subjet kinematics
   // since the symmetry condition depends on the DeltaR between subjets.
-  virtual double symmetry_cut_fn(const PseudoJet & /* p1 */, 
-                                 const PseudoJet & /* p2 */) const;
+  virtual double symmetry_cut_fn(const PseudoJet & p1, 
+                                 const PseudoJet & p2,
+                                 void * optional_R0sqr_ptr = 0) const;
   virtual std::string symmetry_cut_description() const;
 
-private:
+  //private:
   double _beta;         ///< the power of the angular distance to be used
                         ///< in the symmetry condition

external/fastjet/contribs/RecursiveTools/TODO

@@ +1 @@
+For v2.0
+--------
+
+- CHECK: lines 100-102 of RecursiveSymmetryCutBase.hh: this was
+  initially setting the structure members to 0. I think it's highly
+  preferable to set them to their default of -1, signalling the
+  absence of substructure. [0 dR might happen with a perfectly
+  collinear splitting]
+
+- CHECK: are we happy with how the structure info is stored and the
+  recursive_soft_drop_prongs(...) function in RecursiveSoftDrop.hh?
+
+- QUESTION: do we move set_min_deltaR_squared in the base class?
+
+- Add option to define z wrt to original jet?
+
 For future releases?
 --------------------
@@ -7 +23 @@
   auto-indenting. [Consider this for the full contrib??]
 - More generic kinematic cuts?
+- KZ: common base class for IteratedSoftDrop and RecursiveSoftDrop?
+- KZ: make IteratedSoftDrop use Recluster

external/fastjet/contribs/RecursiveTools/VERSION

@@ -1 @@
-1.0.0
+2.0.0-beta2