Index: contrib/contribs/RecursiveTools/tags/2.0.1/TODO =================================================================== --- contrib/contribs/RecursiveTools/tags/2.0.1/TODO (revision 0) +++ contrib/contribs/RecursiveTools/tags/2.0.1/TODO (revision 1283) @@ -0,0 +1,26 @@ +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? +-------------------- +- JDT: In Recluster, change "single" to an enum for user clarity? + [do we have any other option that "sungle" or "subjets"] +- JDT: Should we have a FASTJET_CONTRIB_BEGIN_NAMESPACE? It would + help those of us who use XCode to edit, which is aggressive about + 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 \ No newline at end of file Index: contrib/contribs/RecursiveTools/tags/2.0.1/Recluster.hh =================================================================== --- contrib/contribs/RecursiveTools/tags/2.0.1/Recluster.hh (revision 0) +++ contrib/contribs/RecursiveTools/tags/2.0.1/Recluster.hh (revision 1283) @@ -0,0 +1,183 @@ +#ifndef __FASTJET_CONTRIB_TOOLS_RECLUSTER_HH__ +#define __FASTJET_CONTRIB_TOOLS_RECLUSTER_HH__ + +// $Id$ +// +// Copyright (c) 2014-, Matteo Cacciari, Gavin P. Salam and Gregory Soyez +// +//---------------------------------------------------------------------- +// This file is part of FastJet contrib. +// +// It is free software; you can redistribute it and/or modify it under +// the terms of the GNU General Public License as published by the +// Free Software Foundation; either version 2 of the License, or (at +// your option) any later version. +// +// It is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY +// or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public +// License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this code. If not, see . +//---------------------------------------------------------------------- + +#include +#include // to derive the ReclusterStructure from CompositeJetStructure +#include // to derive Recluster from Transformer +#include +#include + +FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh + +namespace contrib{ + +//---------------------------------------------------------------------- +/// \class Recluster +/// Class that helps reclustering a jet with a new jet definition +/// +/// The result of the reclustering is returned as a single PseudoJet +/// with a CompositeJet structure. The pieces of that PseudoJet will +/// be the individual subjets +/// +/// When constructed from a JetDefinition, that definition will be +/// used to obtain the subjets. When constructed from a JetAlgorithm +/// and parameters (0 parameters for e+e-, just R or R and an extra +/// parameter for others) the recombination scheme will be taken as +/// the same one used to initially cluster the original jet. +/// +/// The result of this transformer depends on its usage. There are two +/// typical use-cases: either we recluster one fat jet into subjets, +/// OR, we recluster the jet with a different jet alg. When Recluster +/// is created from a full jet definition. The last parameter of the +/// constructors below dicatate that behaviour: if "single" is true +/// (the default), a single jet, issued from a regular clustering is +/// returned (if there are more than one, the hardest is taken); +/// otherwise (single==false), the result will be a composite jet with +/// each subjet as pieces +/// +/// Open points for discussion: +/// +/// - do we add an option to force area support? [could be useful +/// e.g. for the filter with a subtractor where area support is +/// mandatory] +/// +class Recluster : public Transformer { +public: + /// define a recluster that decomposes a jet into subjets using a + /// generic JetDefinition + /// + /// \param subjet_def the jet definition applied to obtain the subjets + /// \param single when true, cluster the jet in a single jet (the + /// hardest one) with an associated ClusterSequence, + /// otherwise return a composite jet with subjets + /// as pieces. + Recluster(const JetDefinition & subjet_def, bool single=true) + : _subjet_def(subjet_def), _use_full_def(true), _single(single) {} + + /// define a recluster that decomposes a jet into subjets using a + /// JetAlgorithm and its parameters + /// + /// \param subjet_alg the jet algorithm applied to obtain the subjets + /// \param subjet_radius the jet radius if required + /// \param subjet_extra optional extra parameters for the jet algorithm (only when needed) + /// \param single when true, cluster the jet in a single jet (the + /// hardest one) with an associated ClusterSequence, + /// otherwise return a composite jet with subjets + /// as pieces. + /// + /// Typically, for e+e- algoriothm you should use the third version + /// below with no parameters, for "standard" pp algorithms, just the + /// clustering radius has to be specified and for genkt-type of + /// algorithms, both the radius and the extra parameter have to be + /// specified. + Recluster(JetAlgorithm subjet_alg, double subjet_radius, double subjet_extra, + bool single=true) + : _subjet_alg(subjet_alg), _use_full_def(false), + _subjet_radius(subjet_radius), _has_subjet_radius(true), + _subjet_extra(subjet_extra), _has_subjet_extra(true), _single(single) {} + Recluster(JetAlgorithm subjet_alg, double subjet_radius, bool single=true) + : _subjet_alg(subjet_alg), _use_full_def(false), + _subjet_radius(subjet_radius), _has_subjet_radius(true), + _has_subjet_extra(false), _single(single) {} + Recluster(JetAlgorithm subjet_alg, bool single=true) + : _subjet_alg(subjet_alg), _use_full_def(false), + _has_subjet_radius(false), _has_subjet_extra(false), _single(single) {} + + /// default dtor + virtual ~Recluster(){}; + + //---------------------------------------------------------------------- + // standard Transformer behaviour inherited from the base class + // (i.e. result(), description() and structural info) + + /// runs the reclustering and sets kept and rejected to be the jets of interest + /// (with non-zero rho, they will have been subtracted). + /// + /// \param jet the jet that gets reclustered + /// \return the reclustered jet + virtual PseudoJet result(const PseudoJet & jet) const; + + /// class description + virtual std::string description() const; + + // the type of the associated structure + typedef CompositeJetStructure StructureType; + +private: + /// set the reclustered elements in the simple case of C/A+C/A + void _recluster_cafilt(const std::vector & all_pieces, + std::vector & subjets, + double Rfilt) const; + + /// set the reclustered elements in the generic re-clustering case + void _recluster_generic(const PseudoJet & jet, + std::vector & subjets, + const JetDefinition & subjet_def, + bool do_areas) const; + + // a series of checks + //-------------------------------------------------------------------- + /// get the pieces down to the fundamental pieces + bool _get_all_pieces(const PseudoJet &jet, std::vector &all_pieces) const; + + /// get the common recombiner to all pieces (NULL if none) + const JetDefinition::Recombiner* _get_common_recombiner(const std::vector &all_pieces) const; + + /// construct the proper jet definition ensuring that the recombiner + /// is taken from the underlying pieces (an error is thrown if the + /// pieces do no share a common recombiner) + void _build_jet_def_with_recombiner(const std::vector &all_pieces, + JetDefinition &subjet_def) const; + + /// check if one can apply the simplified trick for C/A subjets + bool _check_ca(const std::vector &all_pieces, + const JetDefinition &subjet_def) const; + + /// check if the jet (or all its pieces) have explicit ghosts + /// (assuming the jet has area support + /// + /// Note that if the jet has an associated cluster sequence that is no + /// longer valid, an error will be thrown + bool _check_explicit_ghosts(const std::vector &all_pieces) const; + + JetDefinition _subjet_def; ///< the jet definition to use to extract the subjets + JetAlgorithm _subjet_alg; ///< the jet algorithm to be used + bool _use_full_def; ///< true when the full JetDefinition is supplied to the ctor + double _subjet_radius; ///< the jet radius (only if needed for the jet alg) + bool _has_subjet_radius; ///< the subjet radius has been specified + double _subjet_extra; ///< the jet alg extra param (only if needed) + bool _has_subjet_extra; ///< the extra param has been specified + + bool _single; ///< (true) return a single jet with a + ///< regular clustering or (false) a + ///< composite jet with subjets as pieces + + static LimitedWarning _explicit_ghost_warning; +}; + +} // namespace contrib + +FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh + +#endif // __FASTJET_CONTRIB_TOOLS_RECLUSTER_HH__ Property changes on: contrib/contribs/RecursiveTools/tags/2.0.1/Recluster.hh ___________________________________________________________________ Added: svn:keywords ## -0,0 +1 ## +Id \ No newline at end of property Index: contrib/contribs/RecursiveTools/tags/2.0.1/IteratedSoftDrop.cc =================================================================== --- contrib/contribs/RecursiveTools/tags/2.0.1/IteratedSoftDrop.cc (revision 0) +++ contrib/contribs/RecursiveTools/tags/2.0.1/IteratedSoftDrop.cc (revision 1283) @@ -0,0 +1,128 @@ +// $Id$ +// +// Copyright (c) 2017-, Jesse Thaler, Kevin Zhou, Gavin P. Salam +// andGregory Soyez +// +// based on arXiv:1704.06266 by Christopher Frye, Andrew J. Larkoski, +// Jesse Thaler, Kevin Zhou +// +//---------------------------------------------------------------------- +// This file is part of FastJet contrib. +// +// It is free software; you can redistribute it and/or modify it under +// the terms of the GNU General Public License as published by the +// Free Software Foundation; either version 2 of the License, or (at +// your option) any later version. +// +// It is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY +// or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public +// License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this code. If not, see . +//---------------------------------------------------------------------- + +#include "IteratedSoftDrop.hh" +#include +#include +#include + +using namespace std; + +FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh + +namespace contrib{ + +//======================================================================== +// implementation of IteratedSoftDropInfo +//======================================================================== + +// returns the angularity with angular exponent alpha and z +// exponent kappa calculated on the zg's and thetag's found by +// iterated SoftDrop +// +// returns 0 if no substructure was found +double IteratedSoftDropInfo::angularity(double alpha, double kappa) const{ + double sum = 0.0; + for (unsigned int i=0; i< _all_zg_thetag.size(); ++i) + sum += pow(_all_zg_thetag[i].first, kappa) * pow(_all_zg_thetag[i].second, alpha); + return sum; +} + +//======================================================================== +// implementation of IteratedSoftDrop +//======================================================================== + +// Constructor. Takes in the standard Soft Drop parameters, an angular cut \f$\theta_{\rm cut}\f$, +// and a choice of angular and symmetry measure. +// +// - beta the Soft Drop beta parameter +// - symmetry_cut the Soft Drop symmetry cut +// - angular_cut the angular cutoff to halt Iterated Soft Drop +// - R0 the angular distance normalization +IteratedSoftDrop::IteratedSoftDrop(double beta, double symmetry_cut, double angular_cut, double R0, + const FunctionOfPseudoJet * subtractor) : + _rsd(beta, symmetry_cut, -1, R0, subtractor){ + _rsd.set_hardest_branch_only(true); + if (angular_cut>0) + _rsd.set_min_deltaR_squared(angular_cut*angular_cut); +} + + +// Full constructor, which takes the following parameters: +// +// \param beta the value of the beta parameter +// \param symmetry_cut the value of the cut on the symmetry measure +// \param symmetry_measure the choice of measure to use to estimate the symmetry +// \param angular_cut the angular cutoff to halt Iterated Soft Drop +// \param R0 the angular distance normalisation [1 by default] +// \param mu_cut the maximal allowed value of mass drop variable mu = m_heavy/m_parent +// \param recursion_choice the strategy used to decide which subjet to recurse into +// \param subtractor an optional pointer to a pileup subtractor (ignored if zero) +IteratedSoftDrop::IteratedSoftDrop(double beta, + double symmetry_cut, + RecursiveSoftDrop::SymmetryMeasure symmetry_measure, + double angular_cut, + double R0, + double mu_cut, + RecursiveSoftDrop::RecursionChoice recursion_choice, + const FunctionOfPseudoJet * subtractor) + : _rsd(beta, symmetry_cut, symmetry_measure, -1, R0, mu_cut, recursion_choice, subtractor){ + _rsd.set_hardest_branch_only(true); + if (angular_cut>0) + _rsd.set_min_deltaR_squared(angular_cut*angular_cut); +} + + +// returns vector of ISD symmetry factors and splitting angles +IteratedSoftDropInfo IteratedSoftDrop::result(const PseudoJet& jet) const{ + PseudoJet rsd_jet = _rsd(jet); + if (! rsd_jet.has_structure_of()) + return IteratedSoftDropInfo(); + return IteratedSoftDropInfo(rsd_jet.structure_of().sorted_zg_and_thetag()); +} + + +std::string IteratedSoftDrop::description() const{ + std::ostringstream oss; + oss << "IteratedSoftDrop with beta =" << _rsd.beta() + << ", symmetry_cut=" << _rsd.symmetry_cut() + << ", R0=" << _rsd.R0(); + + if (_rsd.min_deltaR_squared() >= 0){ + oss << " and angular_cut=" << sqrt(_rsd.min_deltaR_squared()); + } else { + oss << " and no angular_cut"; + } + + if (_rsd.subtractor()){ + oss << ", and with internal subtraction using [" << _rsd.subtractor()->description() << "]"; + } + return oss.str(); +} + + +} // namespace contrib + +FASTJET_END_NAMESPACE Property changes on: contrib/contribs/RecursiveTools/tags/2.0.1/IteratedSoftDrop.cc ___________________________________________________________________ Added: svn:keywords ## -0,0 +1 ## +Id \ No newline at end of property Index: contrib/contribs/RecursiveTools/tags/2.0.1/ModifiedMassDropTagger.cc =================================================================== --- contrib/contribs/RecursiveTools/tags/2.0.1/ModifiedMassDropTagger.cc (revision 0) +++ contrib/contribs/RecursiveTools/tags/2.0.1/ModifiedMassDropTagger.cc (revision 1283) @@ -0,0 +1,43 @@ +// $Id$ +// +// Copyright (c) 2014-, Gavin P. Salam +// +//---------------------------------------------------------------------- +// This file is part of FastJet contrib. +// +// It is free software; you can redistribute it and/or modify it under +// the terms of the GNU General Public License as published by the +// Free Software Foundation; either version 2 of the License, or (at +// your option) any later version. +// +// It is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY +// or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public +// License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this code. If not, see . +//---------------------------------------------------------------------- + +#include "ModifiedMassDropTagger.hh" +#include "fastjet/JetDefinition.hh" +#include "fastjet/ClusterSequenceAreaBase.hh" +#include +#include + +using namespace std; + +FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh + +namespace contrib{ + +//---------------------------------------------------------------------- +string ModifiedMassDropTagger::symmetry_cut_description() const { + ostringstream ostr; + ostr << _symmetry_cut << " [ModifiedMassDropTagger]"; + return ostr.str(); +} + +} // namespace contrib + +FASTJET_END_NAMESPACE Property changes on: contrib/contribs/RecursiveTools/tags/2.0.1/ModifiedMassDropTagger.cc ___________________________________________________________________ Added: svn:keywords ## -0,0 +1 ## +Id \ No newline at end of property Index: contrib/contribs/RecursiveTools/tags/2.0.1/COPYING =================================================================== --- contrib/contribs/RecursiveTools/tags/2.0.1/COPYING (revision 0) +++ contrib/contribs/RecursiveTools/tags/2.0.1/COPYING (revision 1283) @@ -0,0 +1,360 @@ +The RecursiveTools contrib to FastJet is released +under the terms of the GNU General Public License v2 (GPLv2). + +A copy of the GPLv2 is to be found at the end of this file. + +While the GPL license grants you considerable freedom, please bear in +mind that this code's use falls under guidelines similar to those that +are standard for Monte Carlo event generators +(http://www.montecarlonet.org/GUIDELINES). In particular, if you use +this code as part of work towards a scientific publication, whether +directly or contained within another program you should include a +citation to + + arXiv:1307.0007 + arXiv:1402.2657 + arXiv:1704.06266 + +====================================================================== +====================================================================== +====================================================================== + GNU GENERAL PUBLIC LICENSE + Version 2, June 1991 + + Copyright (C) 1989, 1991 Free Software Foundation, Inc. + 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + Everyone is permitted to copy and distribute verbatim copies + of this license document, but changing it is not allowed. + + Preamble + + The licenses for most software are designed to take away your +freedom to share and change it. By contrast, the GNU General Public +License is intended to guarantee your freedom to share and change free +software--to make sure the software is free for all its users. This +General Public License applies to most of the Free Software +Foundation's software and to any other program whose authors commit to +using it. (Some other Free Software Foundation software is covered by +the GNU Library General Public License instead.) You can apply it to +your programs, too. + + When we speak of free software, we are referring to freedom, not +price. Our General Public Licenses are designed to make sure that you +have the freedom to distribute copies of free software (and charge for +this service if you wish), that you receive source code or can get it +if you want it, that you can change the software or use pieces of it +in new free programs; and that you know you can do these things. + + To protect your rights, we need to make restrictions that forbid +anyone to deny you these rights or to ask you to surrender the rights. +These restrictions translate to certain responsibilities for you if you +distribute copies of the software, or if you modify it. + + For example, if you distribute copies of such a program, whether +gratis or for a fee, you must give the recipients all the rights that +you have. You must make sure that they, too, receive or can get the +source code. And you must show them these terms so they know their +rights. + + We protect your rights with two steps: (1) copyright the software, and +(2) offer you this license which gives you legal permission to copy, +distribute and/or modify the software. + + Also, for each author's protection and ours, we want to make certain +that everyone understands that there is no warranty for this free +software. If the software is modified by someone else and passed on, we +want its recipients to know that what they have is not the original, so +that any problems introduced by others will not reflect on the original +authors' reputations. + + Finally, any free program is threatened constantly by software +patents. We wish to avoid the danger that redistributors of a free +program will individually obtain patent licenses, in effect making the +program proprietary. To prevent this, we have made it clear that any +patent must be licensed for everyone's free use or not licensed at all. + + The precise terms and conditions for copying, distribution and +modification follow. + + GNU GENERAL PUBLIC LICENSE + TERMS AND CONDITIONS FOR COPYING, DISTRIBUTION AND MODIFICATION + + 0. This License applies to any program or other work which contains +a notice placed by the copyright holder saying it may be distributed +under the terms of this General Public License. The "Program", below, +refers to any such program or work, and a "work based on the Program" +means either the Program or any derivative work under copyright law: +that is to say, a work containing the Program or a portion of it, +either verbatim or with modifications and/or translated into another +language. (Hereinafter, translation is included without limitation in +the term "modification".) Each licensee is addressed as "you". + +Activities other than copying, distribution and modification are not +covered by this License; they are outside its scope. The act of +running the Program is not restricted, and the output from the Program +is covered only if its contents constitute a work based on the +Program (independent of having been made by running the Program). +Whether that is true depends on what the Program does. + + 1. You may copy and distribute verbatim copies of the Program's +source code as you receive it, in any medium, provided that you +conspicuously and appropriately publish on each copy an appropriate +copyright notice and disclaimer of warranty; keep intact all the +notices that refer to this License and to the absence of any warranty; +and give any other recipients of the Program a copy of this License +along with the Program. + +You may charge a fee for the physical act of transferring a copy, and +you may at your option offer warranty protection in exchange for a fee. + + 2. You may modify your copy or copies of the Program or any portion +of it, thus forming a work based on the Program, and copy and +distribute such modifications or work under the terms of Section 1 +above, provided that you also meet all of these conditions: + + a) You must cause the modified files to carry prominent notices + stating that you changed the files and the date of any change. + + b) You must cause any work that you distribute or publish, that in + whole or in part contains or is derived from the Program or any + part thereof, to be licensed as a whole at no charge to all third + parties under the terms of this License. + + c) If the modified program normally reads commands interactively + when run, you must cause it, when started running for such + interactive use in the most ordinary way, to print or display an + announcement including an appropriate copyright notice and a + notice that there is no warranty (or else, saying that you provide + a warranty) and that users may redistribute the program under + these conditions, and telling the user how to view a copy of this + License. (Exception: if the Program itself is interactive but + does not normally print such an announcement, your work based on + the Program is not required to print an announcement.) + +These requirements apply to the modified work as a whole. If +identifiable sections of that work are not derived from the Program, +and can be reasonably considered independent and separate works in +themselves, then this License, and its terms, do not apply to those +sections when you distribute them as separate works. But when you +distribute the same sections as part of a whole which is a work based +on the Program, the distribution of the whole must be on the terms of +this License, whose permissions for other licensees extend to the +entire whole, and thus to each and every part regardless of who wrote it. + +Thus, it is not the intent of this section to claim rights or contest +your rights to work written entirely by you; rather, the intent is to +exercise the right to control the distribution of derivative or +collective works based on the Program. + +In addition, mere aggregation of another work not based on the Program +with the Program (or with a work based on the Program) on a volume of +a storage or distribution medium does not bring the other work under +the scope of this License. + + 3. You may copy and distribute the Program (or a work based on it, +under Section 2) in object code or executable form under the terms of +Sections 1 and 2 above provided that you also do one of the following: + + a) Accompany it with the complete corresponding machine-readable + source code, which must be distributed under the terms of Sections + 1 and 2 above on a medium customarily used for software interchange; or, + + b) Accompany it with a written offer, valid for at least three + years, to give any third party, for a charge no more than your + cost of physically performing source distribution, a complete + machine-readable copy of the corresponding source code, to be + distributed under the terms of Sections 1 and 2 above on a medium + customarily used for software interchange; or, + + c) Accompany it with the information you received as to the offer + to distribute corresponding source code. (This alternative is + allowed only for noncommercial distribution and only if you + received the program in object code or executable form with such + an offer, in accord with Subsection b above.) + +The source code for a work means the preferred form of the work for +making modifications to it. For an executable work, complete source +code means all the source code for all modules it contains, plus any +associated interface definition files, plus the scripts used to +control compilation and installation of the executable. However, as a +special exception, the source code distributed need not include +anything that is normally distributed (in either source or binary +form) with the major components (compiler, kernel, and so on) of the +operating system on which the executable runs, unless that component +itself accompanies the executable. + +If distribution of executable or object code is made by offering +access to copy from a designated place, then offering equivalent +access to copy the source code from the same place counts as +distribution of the source code, even though third parties are not +compelled to copy the source along with the object code. + + 4. You may not copy, modify, sublicense, or distribute the Program +except as expressly provided under this License. Any attempt +otherwise to copy, modify, sublicense or distribute the Program is +void, and will automatically terminate your rights under this License. +However, parties who have received copies, or rights, from you under +this License will not have their licenses terminated so long as such +parties remain in full compliance. + + 5. You are not required to accept this License, since you have not +signed it. However, nothing else grants you permission to modify or +distribute the Program or its derivative works. These actions are +prohibited by law if you do not accept this License. Therefore, by +modifying or distributing the Program (or any work based on the +Program), you indicate your acceptance of this License to do so, and +all its terms and conditions for copying, distributing or modifying +the Program or works based on it. + + 6. Each time you redistribute the Program (or any work based on the +Program), the recipient automatically receives a license from the +original licensor to copy, distribute or modify the Program subject to +these terms and conditions. You may not impose any further +restrictions on the recipients' exercise of the rights granted herein. +You are not responsible for enforcing compliance by third parties to +this License. + + 7. If, as a consequence of a court judgment or allegation of patent +infringement or for any other reason (not limited to patent issues), +conditions are imposed on you (whether by court order, agreement or +otherwise) that contradict the conditions of this License, they do not +excuse you from the conditions of this License. If you cannot +distribute so as to satisfy simultaneously your obligations under this +License and any other pertinent obligations, then as a consequence you +may not distribute the Program at all. For example, if a patent +license would not permit royalty-free redistribution of the Program by +all those who receive copies directly or indirectly through you, then +the only way you could satisfy both it and this License would be to +refrain entirely from distribution of the Program. + +If any portion of this section is held invalid or unenforceable under +any particular circumstance, the balance of the section is intended to +apply and the section as a whole is intended to apply in other +circumstances. + +It is not the purpose of this section to induce you to infringe any +patents or other property right claims or to contest validity of any +such claims; this section has the sole purpose of protecting the +integrity of the free software distribution system, which is +implemented by public license practices. Many people have made +generous contributions to the wide range of software distributed +through that system in reliance on consistent application of that +system; it is up to the author/donor to decide if he or she is willing +to distribute software through any other system and a licensee cannot +impose that choice. + +This section is intended to make thoroughly clear what is believed to +be a consequence of the rest of this License. + + 8. If the distribution and/or use of the Program is restricted in +certain countries either by patents or by copyrighted interfaces, the +original copyright holder who places the Program under this License +may add an explicit geographical distribution limitation excluding +those countries, so that distribution is permitted only in or among +countries not thus excluded. In such case, this License incorporates +the limitation as if written in the body of this License. + + 9. The Free Software Foundation may publish revised and/or new versions +of the General Public License from time to time. Such new versions will +be similar in spirit to the present version, but may differ in detail to +address new problems or concerns. + +Each version is given a distinguishing version number. If the Program +specifies a version number of this License which applies to it and "any +later version", you have the option of following the terms and conditions +either of that version or of any later version published by the Free +Software Foundation. If the Program does not specify a version number of +this License, you may choose any version ever published by the Free Software +Foundation. + + 10. If you wish to incorporate parts of the Program into other free +programs whose distribution conditions are different, write to the author +to ask for permission. For software which is copyrighted by the Free +Software Foundation, write to the Free Software Foundation; we sometimes +make exceptions for this. Our decision will be guided by the two goals +of preserving the free status of all derivatives of our free software and +of promoting the sharing and reuse of software generally. + + NO WARRANTY + + 11. BECAUSE THE PROGRAM IS LICENSED FREE OF CHARGE, THERE IS NO WARRANTY +FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE LAW. EXCEPT WHEN +OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR OTHER PARTIES +PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESSED +OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF +MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS +TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE +PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY SERVICING, +REPAIR OR CORRECTION. + + 12. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING +WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY MODIFY AND/OR +REDISTRIBUTE THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, +INCLUDING ANY GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING +OUT OF THE USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED +TO LOSS OF DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY +YOU OR THIRD PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER +PROGRAMS), EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE +POSSIBILITY OF SUCH DAMAGES. + + END OF TERMS AND CONDITIONS + + How to Apply These Terms to Your New Programs + + If you develop a new program, and you want it to be of the greatest +possible use to the public, the best way to achieve this is to make it +free software which everyone can redistribute and change under these terms. + + To do so, attach the following notices to the program. It is safest +to attach them to the start of each source file to most effectively +convey the exclusion of warranty; and each file should have at least +the "copyright" line and a pointer to where the full notice is found. + + + Copyright (C) + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program; if not, write to the Free Software + Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + + +Also add information on how to contact you by electronic and paper mail. + +If the program is interactive, make it output a short notice like this +when it starts in an interactive mode: + + Gnomovision version 69, Copyright (C) year name of author + Gnomovision comes with ABSOLUTELY NO WARRANTY; for details type `show w'. + This is free software, and you are welcome to redistribute it + under certain conditions; type `show c' for details. + +The hypothetical commands `show w' and `show c' should show the appropriate +parts of the General Public License. Of course, the commands you use may +be called something other than `show w' and `show c'; they could even be +mouse-clicks or menu items--whatever suits your program. + +You should also get your employer (if you work as a programmer) or your +school, if any, to sign a "copyright disclaimer" for the program, if +necessary. Here is a sample; alter the names: + + Yoyodyne, Inc., hereby disclaims all copyright interest in the program + `Gnomovision' (which makes passes at compilers) written by James Hacker. + + , 1 April 1989 + Ty Coon, President of Vice + +This General Public License does not permit incorporating your program into +proprietary programs. If your program is a subroutine library, you may +consider it more useful to permit linking proprietary applications with the +library. If this is what you want to do, use the GNU Library General +Public License instead of this License. Index: contrib/contribs/RecursiveTools/tags/2.0.1/example_mmdt_ee.cc =================================================================== --- contrib/contribs/RecursiveTools/tags/2.0.1/example_mmdt_ee.cc (revision 0) +++ contrib/contribs/RecursiveTools/tags/2.0.1/example_mmdt_ee.cc (revision 1283) @@ -0,0 +1,142 @@ +//---------------------------------------------------------------------- +/// \file example_ee.cc +/// +/// This example program is meant to illustrate how to use +/// RecursiveTools for e+e- events. It is done using the +/// ModifiedMassDropTagger class but the same strategy would work as +/// well for SoftDrop, RecursiveSoftDrop and IteratedSoftDrop +/// +/// Run this example with +/// +/// \verbatim +/// ./example_ee < ../data/single-ee-event.dat +/// \endverbatim +//---------------------------------------------------------------------- + +// $Id$ +// +// Copyright (c) 2017-, Gavin P. Salam, Gregory Soyez, Jesse Thaler, +// Kevin Zhou, Frederic Dreyer +// +//---------------------------------------------------------------------- +// This file is part of FastJet contrib. +// +// It is free software; you can redistribute it and/or modify it under +// the terms of the GNU General Public License as published by the +// Free Software Foundation; either version 2 of the License, or (at +// your option) any later version. +// +// It is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY +// or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public +// License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this code. If not, see . +//---------------------------------------------------------------------- + +#include +#include + +#include +#include +#include +#include "fastjet/ClusterSequence.hh" +#include "fastjet/tools/Filter.hh" +#include "ModifiedMassDropTagger.hh" // In external code, this should be fastjet/contrib/ModifiedMassDropTagger.hh + +using namespace std; +using namespace fastjet; + +// forward declaration to make things clearer +void read_event(vector &event); +ostream & operator<<(ostream &, const PseudoJet &); + +//---------------------------------------------------------------------- +int main(){ + + //---------------------------------------------------------- + // read in input particles + vector event; + read_event(event); + cout << "# read an event with " << event.size() << " particles" << endl; + + // first get some Cambridge/Aachen jets + double R = 1.0; + JetDefinition jet_def(ee_genkt_algorithm, R, 0.0); + ClusterSequence cs(event, jet_def); + + double Emin = 10.0; + Selector sel_jets = SelectorEMin(Emin); + vector jets = sorted_by_E(sel_jets(cs.inclusive_jets())); + + // give the tagger a short name + typedef contrib::ModifiedMassDropTagger MMDT; + + // This version uses the following setup: + // - use energy for the symmetry measure + // Note: since the mMDT does not require angular information, + // here we could use either theta_E or cos_theta_E (it would + // only change the value of DeltaR). For + // SoftDrop/RecursiveSoftDrop/IteratedSoftDrop, this would make + // a difference and we have + // DeltaR_{ij}^2 = theta_{ij}^2 (theta_E) + // DeltaR_{ij}^2 = 2 [1-cos(theta_{ij}^2)] (cos_theta_E) + // + // - recurse into the branch with the largest energy + // + // - use a symmetry cut, with no mass-drop requirement + double z_cut = 0.20; + MMDT tagger(z_cut, + MMDT::cos_theta_E, + std::numeric_limits::infinity(), + MMDT::larger_E); + cout << "tagger is: " << tagger.description() << endl; + + for (unsigned ijet = 0; ijet < jets.size(); ijet++) { + // first run MMDT and examine the output + PseudoJet tagged_jet = tagger(jets[ijet]); + cout << endl; + cout << "original jet: " << jets[ijet] << endl; + cout << "tagged jet: " << tagged_jet << endl; + if (tagged_jet == 0) continue; // If symmetry condition not satisified, jet is not tagged + cout << " delta_R between subjets: " << tagged_jet.structure_of().delta_R() << endl; + cout << " symmetry measure(z): " << tagged_jet.structure_of().symmetry() << endl; + cout << " mass drop(mu): " << tagged_jet.structure_of().mu() << endl; + + cout << endl; + } + + return 0; +} + +//---------------------------------------------------------------------- +/// read in input particles +void read_event(vector &event){ + string line; + while (getline(cin, line)) { + istringstream linestream(line); + // take substrings to avoid problems when there are extra "pollution" + // characters (e.g. line-feed). + if (line.substr(0,4) == "#END") {return;} + if (line.substr(0,1) == "#") {continue;} + double px,py,pz,E; + linestream >> px >> py >> pz >> E; + PseudoJet particle(px,py,pz,E); + + // push event onto back of full_event vector + event.push_back(particle); + } +} + +//---------------------------------------------------------------------- +/// overloaded jet info output +ostream & operator<<(ostream & ostr, const PseudoJet & jet) { + if (jet == 0) { + ostr << " 0 "; + } else { + ostr << " E = " << jet.pt() + << " m = " << jet.m(); + } + return ostr; +} Property changes on: contrib/contribs/RecursiveTools/tags/2.0.1/example_mmdt_ee.cc ___________________________________________________________________ Added: svn:keywords ## -0,0 +1 ## +Id \ No newline at end of property Index: contrib/contribs/RecursiveTools/tags/2.0.1/RecursiveSoftDrop.cc =================================================================== --- contrib/contribs/RecursiveTools/tags/2.0.1/RecursiveSoftDrop.cc (revision 0) +++ contrib/contribs/RecursiveTools/tags/2.0.1/RecursiveSoftDrop.cc (revision 1283) @@ -0,0 +1,475 @@ +// $Id$ +// +// Copyright (c) 2017-, Gavin P. Salam, Gregory Soyez, Jesse Thaler, +// Kevin Zhou, Frederic Dreyer +// +//---------------------------------------------------------------------- +// This file is part of FastJet contrib. +// +// It is free software; you can redistribute it and/or modify it under +// the terms of the GNU General Public License as published by the +// Free Software Foundation; either version 2 of the License, or (at +// your option) any later version. +// +// It is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY +// or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public +// License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this code. If not, see . +//---------------------------------------------------------------------- + +#include "RecursiveSoftDrop.hh" +#include "fastjet/ClusterSequence.hh" + +using namespace std; + +FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh + +namespace contrib{ + +namespace internal_recursive_softdrop{ + +//======================================================================== +/// \class RSDHistoryElement +/// a helper class to help keeping track od the RSD tree +/// +/// The element is created at the top of a branch and updated each +/// time one grooms something away. +class RSDHistoryElement{ +public: + RSDHistoryElement(const PseudoJet &jet, const RecursiveSoftDrop* rsd_ptr, double R0sqr) : + R0_squared(R0sqr), + child1_in_history(-1), child2_in_history(-1), symmetry(-1.0), mu2(-1.0){ + reset(jet, rsd_ptr); + } + + void reset(const PseudoJet &jet, const RecursiveSoftDrop* rsd_ptr){ + current_in_ca_tree = jet.cluster_hist_index(); + PseudoJet piece1, piece2; + theta_squared = (jet.has_parents(piece1, piece2)) + ? rsd_ptr->squared_geometric_distance(piece1,piece2) : 0.0; + } + + int current_in_ca_tree; ///< (history) index of the current particle in the C/A tree + double theta_squared; ///< squared angle at which this decays + double R0_squared; ///< squared angle at the top of the branch + ///< (used for RSD with dynamic_R0) + int child1_in_history; ///< hardest of the 2 decay products (-1 if untagged) + int child2_in_history; ///< softest of the 2 decay products (-1 if untagged) + + // info about what has been dropped and the local substructure + vector dropped_delta_R; + vector dropped_symmetry; + vector dropped_mu; + double symmetry, mu2; +}; + + +/// \class OrderRSDHistoryElements +/// angular ordering of (pointers to) the history elements +/// +/// our priority queue will use pointers to these elements that are +/// ordered in angle (of the objects they point to) +class OrderRSDHistoryElements{ +public: + bool operator()(const RSDHistoryElement *e1, const RSDHistoryElement *e2) const { + return e1->theta_squared < e2->theta_squared; + } +}; + +} // internal_recursive_softdrop + +//======================================================================== + +// initialise all the flags and parameters to their default value +void RecursiveSoftDrop::set_defaults(){ + set_fixed_depth_mode(false); + set_dynamical_R0(false); + set_hardest_branch_only(false); + set_min_deltaR_squared(-1.0); +} + +// description of the tool +string RecursiveSoftDrop::description() const{ + ostringstream res; + res << "recursive application of [" + << RecursiveSymmetryCutBase::description() + << "]"; + + if (_fixed_depth){ + res << ", recursively applied down to a maximal depth of N="; + if (_n==-1) res << "infinity"; else res << _n; + } else { + res << ", applied N="; + if (_n==-1) res << "infinity"; else res << _n; + res << " times"; + } + + if (_dynamical_R0) + res << ", with R0 dynamically scaled"; + else + res << ", with R0 kept fixed"; + + if (_hardest_branch_only) + res << ", following only the hardest branch"; + + if (_min_dR2>0) + res << ", with minimal angle (squared) = " << _min_dR2; + + return res.str(); +} + + +// action on a single jet with RecursiveSoftDrop. +// +// uses "result_fixed_tags" by default (i.e. recurse from R0 to +// smaller angles until n SD conditions have been met), or +// "result_fixed_depth" where each of the previous SD branches are +// recirsed into down to a depth of n. +PseudoJet RecursiveSoftDrop::result(const PseudoJet &jet) const{ + return _fixed_depth ? result_fixed_depth(jet) : result_fixed_tags(jet); +} + +// this routine applies the Soft Drop criterion recursively on the +// CA tree until we find n subjets (or until it converges), and +// adds them together into a groomed PseudoJet +PseudoJet RecursiveSoftDrop::result_fixed_tags(const PseudoJet &jet) const { + // start by reclustering jet with C/A algorithm + PseudoJet ca_jet = _recluster_if_needed(jet); + + if (! ca_jet.has_valid_cluster_sequence()){ + throw Error("RecursiveSoftDrop can only be applied to jets associated to a (valid) cluster sequence"); + } + + const ClusterSequence *cs = ca_jet.validated_cluster_sequence(); + const vector &cs_history = cs->history(); + const vector &cs_jets = cs->jets(); + + // initialise counter to 1 subjet (i.e. the full ca_jet) + int n_tagged = 0; + int max_njet = ca_jet.constituents().size(); + + // create the list of branches + unsigned int max_history_size = 2*max_njet; + if ((_n>0) && (_n history; + history.reserve(max_history_size); // could be one shorter + history.push_back(internal_recursive_softdrop::RSDHistoryElement(ca_jet, this, _R0sqr)); + + // create a priority queue containing the subjets and a comparison definition + priority_queue, internal_recursive_softdrop::OrderRSDHistoryElements> active_branches; + active_branches.push(& (history[0])); + + PseudoJet parent, piece1, piece2; + double sym, mu2; + + // loop over C/A tree until we reach the appropriate number of subjets + while ((continue_grooming(n_tagged)) && (active_branches.size())) { + // get the element corresponding to the max dR and the associated PJ + internal_recursive_softdrop::RSDHistoryElement * elm = active_branches.top(); + PseudoJet parent = cs_jets[cs_history[elm->current_in_ca_tree].jetp_index]; + + // do one step of SD + RecursionStatus status = recurse_one_step(parent, piece1, piece2, sym, mu2, &elm->R0_squared); + + // check if we passed the SD condition + if (status==recursion_success){ + // check for the optional angular cut + if ((_min_dR2 > 0) && (squared_geometric_distance(piece1,piece2) < _min_dR2)) + break; + + // both subjets are kept in the list for potential further de-clustering + elm->child1_in_history = history.size(); + elm->child2_in_history = history.size()+1; + elm->symmetry = sym; + elm->mu2 = mu2; + active_branches.pop(); + + // update the history + double next_R0_squared = (_dynamical_R0) + ? piece1.squared_distance(piece2) : elm->R0_squared; + + internal_recursive_softdrop::RSDHistoryElement elm1(piece1, this, next_R0_squared); + history.push_back(elm1); + active_branches.push(&(history.back())); + internal_recursive_softdrop::RSDHistoryElement elm2(piece2, this, next_R0_squared); + history.push_back(elm2); + if (!_hardest_branch_only){ + active_branches.push(&(history.back())); + } + + ++n_tagged; + } else if (status==recursion_dropped){ + // check for the optional angular cut + if ((_min_dR2 > 0) && (squared_geometric_distance(piece1,piece2) < _min_dR2)) + break; + + active_branches.pop(); + // tagging failed and the softest branch should be dropped + // keep track of what has been groomed away + max_njet -= piece2.constituents().size(); + elm->dropped_delta_R .push_back((elm->theta_squared >= 0) ? sqrt(elm->theta_squared) : -sqrt(elm->theta_squared)); + elm->dropped_symmetry.push_back(sym); + elm->dropped_mu .push_back((mu2>=0) ? sqrt(mu2) : -sqrt(mu2)); + + // keep the hardest branch in the recursion + elm->reset(piece1, this); + active_branches.push(elm); + } else if (status==recursion_no_parents){ + if (_min_dR2 > 0) break; + active_branches.pop(); + // nothing specific to do: we just keep the curent jet as a "leaf" + } else { // recursion_issue + active_branches.pop(); + // we've met an issue + // if the piece2 is null as well, it means we've had a critical problem. + // In that case, return an empty PseudoJet + if (piece2 == 0) return PseudoJet(); + + // otherwise, we should consider "piece2" as a final particle + // not to be recursed into + if (_min_dR2 > 0) break; + max_njet -= (piece2.constituents().size()-1); + break; + } + + // If the missing number of tags is exactly the number of objects + // we have left in the recursion, stop + if (n_tagged == max_njet) break; + } + + // now we have a bunch of history elements that we can use to build the final jet + vector mapped_to_history(history.size()); + unsigned int history_index = history.size(); + do { + --history_index; + const internal_recursive_softdrop::RSDHistoryElement & elm = history[history_index]; + + // two kinds of events: either just a final leave, potentially with grooming + // or a brandhing (also with potential grooming at the end) + if (elm.child1_in_history<0){ + // this is a leaf, i.e. with no further substructure + PseudoJet & subjet = mapped_to_history[history_index] + = cs_jets[cs_history[elm.current_in_ca_tree].jetp_index]; + + StructureType * structure = new StructureType(subjet); + if (has_verbose_structure()){ + structure->set_verbose(true); + structure->set_dropped_delta_R (elm.dropped_delta_R); + structure->set_dropped_symmetry(elm.dropped_symmetry); + structure->set_dropped_mu (elm.dropped_mu); + } + subjet.set_structure_shared_ptr(SharedPtr(structure)); + } else { + PseudoJet & subjet = mapped_to_history[history_index] + = join(mapped_to_history[elm.child1_in_history], mapped_to_history[elm.child2_in_history]); + StructureType * structure = new StructureType(subjet, sqrt(elm.theta_squared), elm.symmetry, sqrt(elm.mu2)); + if (has_verbose_structure()){ + structure->set_verbose(true); + structure->set_dropped_delta_R (elm.dropped_delta_R); + structure->set_dropped_symmetry(elm.dropped_symmetry); + structure->set_dropped_mu (elm.dropped_mu); + } + subjet.set_structure_shared_ptr(SharedPtr(structure)); + } + } while (history_index>0); + + return mapped_to_history[0]; +} + +// this routine applies the Soft Drop criterion recursively on the +// CA tree, recursing into all the branches found during the previous iteration +// until n layers have been found (or until it converges) +PseudoJet RecursiveSoftDrop::result_fixed_depth(const PseudoJet &jet) const { + // start by reclustering jet with C/A algorithm + PseudoJet ca_jet = _recluster_if_needed(jet); + + if (! ca_jet.has_valid_cluster_sequence()){ + throw Error("RecursiveSoftDrop can only be applied to jets associated to a (valid) cluster sequence"); + } + + const ClusterSequence *cs = ca_jet.validated_cluster_sequence(); + const vector &cs_history = cs->history(); + const vector &cs_jets = cs->jets(); + + // initialise counter to 1 subjet (i.e. the full ca_jet) + int n_depth = 0; + int max_njet = ca_jet.constituents().size(); + + // create the list of branches + unsigned int max_history_size = 2*max_njet; + //if ((_n>0) && (_n history; + history.reserve(max_history_size); // could be one shorter + history.push_back(internal_recursive_softdrop::RSDHistoryElement(ca_jet, this, _R0sqr)); + history.back().theta_squared = _R0sqr; + + // create a priority queue containing the subjets and a comparison definition + list active_branches; + active_branches.push_back(& (history[0])); + + PseudoJet parent, piece1, piece2; + + while ((continue_grooming(n_depth)) && (active_branches.size())) { + // loop over all the branches and look for substructure + list::iterator hist_it=active_branches.begin(); + while (hist_it!=active_branches.end()){ + // get the element corresponding to the max dR and the associated PJ + internal_recursive_softdrop::RSDHistoryElement * elm = (*hist_it); + PseudoJet parent = cs_jets[cs_history[elm->current_in_ca_tree].jetp_index]; + + // we need to iterate this branch until we find some substructure + PseudoJet result_sd; + if (_dynamical_R0){ + SoftDrop sd(_beta, _symmetry_cut, symmetry_measure(), sqrt(elm->theta_squared), + mu_cut(), recursion_choice(), subtractor()); + sd.set_reclustering(false); + sd.set_verbose_structure(has_verbose_structure()); + result_sd = sd(parent); + } else { + result_sd = SoftDrop::result(parent); + } + + // if we had an empty PJ, that means we ran into some problems. + // just return an empty PJ ourselves + if (result_sd == 0) return PseudoJet(); + + // update the history element to reflect our iteration + elm->current_in_ca_tree = result_sd.cluster_hist_index(); + + if (has_verbose_structure()){ + elm->dropped_delta_R = result_sd.structure_of().dropped_delta_R(); + elm->dropped_symmetry = result_sd.structure_of().dropped_symmetry(); + elm->dropped_mu = result_sd.structure_of().dropped_mu(); + } + + // if some substructure was found: + if (result_sd.structure_of().has_substructure()){ + // update the history element to reflect our iteration + elm->child1_in_history = history.size(); + elm->child2_in_history = history.size()+1; + elm->theta_squared = result_sd.structure_of().delta_R(); + elm->theta_squared *= elm->theta_squared; + elm->symmetry = result_sd.structure_of().symmetry(); + elm->mu2 = result_sd.structure_of().mu(); + elm->mu2 *= elm->mu2; + + // the next iteration will have to handle 2 new history + // elements (the R0squared argument here is unused) + result_sd.has_parents(piece1, piece2); + internal_recursive_softdrop::RSDHistoryElement elm1(piece1, this, _R0sqr); + history.push_back(elm1); + // insert it in the active branches if needed + if (elm1.theta_squared>0) + active_branches.insert(hist_it,&(history.back())); // insert just before + + internal_recursive_softdrop::RSDHistoryElement elm2(piece2, this, _R0sqr); + history.push_back(elm2); + if ((!_hardest_branch_only) && (elm2.theta_squared>0)){ + active_branches.insert(hist_it,&(history.back())); // insert just before + } + } + // otherwise we've just reached the end of the recursion the + // history information has been updated above + // + // we just need to make sure that we do not recurse into that + // element any longer + + list::iterator current = hist_it; + ++hist_it; + active_branches.erase(current); + } // loop over branches at current depth + ++n_depth; + } // loop over depth + + // now we have a bunch of history elements that we can use to build the final jet + vector mapped_to_history(history.size()); + unsigned int history_index = history.size(); + do { + --history_index; + const internal_recursive_softdrop::RSDHistoryElement & elm = history[history_index]; + + // two kinds of events: either just a final leave, poteitially with grooming + // or a brandhing (also with potential grooming at the end) + if (elm.child1_in_history<0){ + // this is a leaf, i.e. with no further sustructure + PseudoJet & subjet = mapped_to_history[history_index] + = cs_jets[cs_history[elm.current_in_ca_tree].jetp_index]; + + StructureType * structure = new StructureType(subjet); + if (has_verbose_structure()){ + structure->set_verbose(true); + } + subjet.set_structure_shared_ptr(SharedPtr(structure)); + } else { + PseudoJet & subjet = mapped_to_history[history_index] + = join(mapped_to_history[elm.child1_in_history], mapped_to_history[elm.child2_in_history]); + StructureType * structure = new StructureType(subjet, sqrt(elm.theta_squared), elm.symmetry, sqrt(elm.mu2)); + if (has_verbose_structure()){ + structure->set_verbose(true); + structure->set_dropped_delta_R (elm.dropped_delta_R); + structure->set_dropped_symmetry(elm.dropped_symmetry); + structure->set_dropped_mu (elm.dropped_mu); + } + subjet.set_structure_shared_ptr(SharedPtr(structure)); + } + } while (history_index>0); + + return mapped_to_history[0]; +} + + +//======================================================================== +// implementation of the helpers +//======================================================================== + +// helper to get all the prongs in a jet that has been obtained using +// RecursiveSoftDrop (instead of recursively parsing the 1->2 +// composite jet structure) +vector recursive_soft_drop_prongs(const PseudoJet & rsd_jet){ + // make sure that the jet has the appropriate RecursiveSoftDrop structure + if (!rsd_jet.has_structure_of()) + return vector(); + + // if this jet has no substructure, just return a 1-prong object + if (!rsd_jet.structure_of().has_substructure()) + return vector(1, rsd_jet); + + // otherwise fill a vector with all the prongs (no specific ordering) + vector prongs; + + // parse the list of PseudoJet we still need to deal with + vector to_parse = rsd_jet.pieces(); // valid both for a C/A recombination step or a RSD join + unsigned int i_parse = 0; + while (i_parse()) && + (current.structure_of().has_substructure())){ + // if this has some deeper substructure, add it to the list of + // things to further process + vector pieces = current.pieces(); + assert(pieces.size() == 2); + to_parse[i_parse] = pieces[0]; + to_parse.push_back(pieces[1]); + } else { + // no further substructure, just add this as a branch + prongs.push_back(current); + ++i_parse; + } + } + + return prongs; +} + +} + +FASTJET_END_NAMESPACE Property changes on: contrib/contribs/RecursiveTools/tags/2.0.1/RecursiveSoftDrop.cc ___________________________________________________________________ Added: svn:keywords ## -0,0 +1 ## +Id \ No newline at end of property Index: contrib/contribs/RecursiveTools/tags/2.0.1/example_recursive_softdrop.ref =================================================================== --- contrib/contribs/RecursiveTools/tags/2.0.1/example_recursive_softdrop.ref (revision 0) +++ contrib/contribs/RecursiveTools/tags/2.0.1/example_recursive_softdrop.ref (revision 1283) @@ -0,0 +1,83 @@ +# read an event with 354 particles +#-------------------------------------------------------------------------- +# FastJet release 3.3.1-devel +# M. Cacciari, G.P. Salam and G. Soyez +# A software package for jet finding and analysis at colliders +# http://fastjet.fr +# +# Please cite EPJC72(2012)1896 [arXiv:1111.6097] if you use this package +# for scientific work and optionally PLB641(2006)57 [hep-ph/0512210]. +# +# FastJet is provided without warranty under the terms of the GNU GPLv2. +# It uses T. Chan's closest pair algorithm, S. Fortune's Voronoi code +# and 3rd party plugin jet algorithms. See COPYING file for details. +#-------------------------------------------------------------------------- +RecursiveSoftDrop groomer is: recursive application of [Recursive Groomer with a symmetry cut scalar_z > 0.2 (theta/1)^0.5 [SoftDrop], no mass-drop requirement, recursion into the subjet with larger pt], applied N=4 times, with R0 dynamically scaled + +original jet: pt = 983.387 m = 39.9912 y = -0.867307 phi = 2.90511 +RecursiveSoftDropped jet: pt = 811.261 m = 6.45947 y = -0.87094 phi = 2.9083 + +Prongs with clustering information +---------------------------------- + branch branch N_groomed max loc substructure + pt mass loc tot zdrop zg thetag + +--> 811.261 6.45947 13 13 0.0294298 0.154333 0.0200353 + +--> 669.635 1.98393 2 2 0 + +--> 141.627 0.765329 0 0 0 0.199906 0.00734861 + +--> 113.316 0.462946 0 0 0 0.208947 0.00550512 + | +--> 89.6387 0.211913 0 0 0 + | +--> 23.6769 0.13957 0 0 0 + +--> 28.3123 0.170026 0 0 0 0.248173 0.00447786 + +--> 21.286 0.13957 0 0 0 + +--> 7.02636 -3.21461e-05 0 0 0 + +Prongs without clustering information +------------------------------------- +(Raw) list of prongs: + pt mass + 0 669.635 1.98393 + 1 89.6387 0.211913 + 2 21.286 0.13957 + 3 23.6769 0.13957 + 4 7.02636 -3.21461e-05 + +Groomed prongs information: +index zg thetag + 1 0.154333 0.0200353 + 2 0.199906 0.00734861 + 3 0.208947 0.00550512 + 4 0.248173 0.00447786 + +original jet: pt = 908.098 m = 87.7124 y = 0.219482 phi = 6.03487 +RecursiveSoftDropped jet: pt = 830.517 m = 4.91035 y = 0.223054 phi = 6.02995 + +Prongs with clustering information +---------------------------------- + branch branch N_groomed max loc substructure + pt mass loc tot zdrop zg thetag + +--> 830.517 4.91035 12 13 0.0232056 0.060784 0.0153863 + +--> 778.731 3.66481 0 1 0 0.235041 0.0101382 + | +--> 599.106 0.403809 1 1 0 + | +--> 179.628 0.853363 0 1 0 0.25773 0.00871739 + | +--> 131.15 0.378456 1 1 0.0606587 0.315246 0.00417298 + | | +--> 89.8058 0.107191 0 0 0 + | | +--> 41.3448 0.13957 0 0 0 + | +--> 48.4785 0.13957 0 0 0 + +--> 51.7916 0.13957 0 0 0 + +Prongs without clustering information +------------------------------------- +(Raw) list of prongs: + pt mass + 0 599.106 0.403809 + 1 51.7916 0.13957 + 2 89.8058 0.107191 + 3 48.4785 0.13957 + 4 41.3448 0.13957 + +Groomed prongs information: +index zg thetag + 1 0.060784 0.0153863 + 2 0.235041 0.0101382 + 3 0.25773 0.00871739 + 4 0.315246 0.00417298 Index: contrib/contribs/RecursiveTools/tags/2.0.1/BottomUpSoftDrop.cc =================================================================== --- contrib/contribs/RecursiveTools/tags/2.0.1/BottomUpSoftDrop.cc (revision 0) +++ contrib/contribs/RecursiveTools/tags/2.0.1/BottomUpSoftDrop.cc (revision 1283) @@ -0,0 +1,324 @@ +// $Id$ +// +// Copyright (c) 2017-, Gavin P. Salam, Gregory Soyez, Jesse Thaler, +// Kevin Zhou, Frederic Dreyer +// +//---------------------------------------------------------------------- +// This file is part of FastJet contrib. +// +// It is free software; you can redistribute it and/or modify it under +// the terms of the GNU General Public License as published by the +// Free Software Foundation; either version 2 of the License, or (at +// your option) any later version. +// +// It is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY +// or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public +// License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this code. If not, see . +//---------------------------------------------------------------------- + +#include "BottomUpSoftDrop.hh" +#include +#include +#include +#include +#include "fastjet/ClusterSequenceActiveAreaExplicitGhosts.hh" +#include "fastjet/Selector.hh" +#include "fastjet/config.h" + + +using namespace std; + + +FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh + +namespace contrib{ + +//---------------------------------------------------------------------- +// BottomUpSoftDrop class +//---------------------------------------------------------------------- + +// action on a single jet +PseudoJet BottomUpSoftDrop::result(const PseudoJet &jet) const{ + // soft drop can only be applied to jets that have constituents + if (!jet.has_constituents()){ + throw Error("BottomUpSoftDrop: trying to apply the Soft Drop transformer to a jet that has no constituents"); + } + + // if the jet has area support and there are explicit ghosts, we can + // transfer that support to the internal re-clustering + // + // Note that this is just meant to maintain the information since + // all the jes will have a 0 area + bool do_areas = jet.has_area() && _check_explicit_ghosts(jet); + + // build the soft drop plugin + BottomUpSoftDropPlugin * softdrop_plugin; + + // for some constructors, we get the recombiner from the + // input jet -- some acrobatics are needed + if (_get_recombiner_from_jet) { + JetDefinition jet_def = _jet_def; + + // if all the pieces have a shared recombiner, we'll use that + // one. Otherwise, use the one from _jet_def as a fallback. + JetDefinition jet_def_for_recombiner; + if (_check_common_recombiner(jet, jet_def_for_recombiner)){ +#if FASTJET_VERSION_NUMBER >= 30100 + // Note that this is better than the option directly passing the + // recombiner (for cases where th ejet def own its recombiner) + // but it's only available in FJ>=3.1 + jet_def.set_recombiner(jet_def_for_recombiner); +#else + jet_def.set_recombiner(jet_def_for_recombiner.recombiner()); +#endif + } + softdrop_plugin = new BottomUpSoftDropPlugin(jet_def, _beta, _symmetry_cut, _R0); + } else { + softdrop_plugin = new BottomUpSoftDropPlugin(_jet_def, _beta, _symmetry_cut, _R0); + } + + // now recluster the constituents of the jet with that plugin + JetDefinition internal_jet_def(softdrop_plugin); + // flag the plugin for automatic deletion _before_ we make + // copies (so that as long as the copies are also present + // it doesn't get deleted). + internal_jet_def.delete_plugin_when_unused(); + + ClusterSequence * cs; + if (do_areas){ + vector particles, ghosts; + SelectorIsPureGhost().sift(jet.constituents(), ghosts, particles); + // determine the ghost area from the 1st ghost (if none, any value + // will do, as the area will be 0 and subtraction will have + // no effect!) + double ghost_area = (ghosts.size()) ? ghosts[0].area() : 0.01; + cs = new ClusterSequenceActiveAreaExplicitGhosts(particles, internal_jet_def, + ghosts, ghost_area); + } else { + cs = new ClusterSequence(jet.constituents(), internal_jet_def); + } + + PseudoJet result_local = SelectorNHardest(1)(cs->inclusive_jets())[0]; + BottomUpSoftDropStructure * s = new BottomUpSoftDropStructure(result_local); + s->_beta = _beta; + s->_symmetry_cut = _symmetry_cut; + s->_R0 = _R0; + result_local.set_structure_shared_ptr(SharedPtr(s)); + + // make sure things remain persistent -- i.e. tell the jet definition + // and the cluster sequence that it is their responsibility to clean + // up memory once the "result" reaches the end of its life in the user's + // code. (The CS deletes itself when the result goes out of scope and + // that also triggers deletion of the plugin) + cs->delete_self_when_unused(); + + return result_local; +} + +// global grooming on a full event +// note: does not support jet areas +vector BottomUpSoftDrop::global_grooming(const vector & event) const { + // start by reclustering the event into one very large jet + ClusterSequence cs(event, _jet_def); + std::vector global_jet = SelectorNHardest(1)(cs.inclusive_jets()); + // if the event is empty, do nothing + if (global_jet.size() == 0) return vector(); + + // apply the groomer to the large jet + PseudoJet result = this->result(global_jet[0]); + return result.constituents(); +} + +// check if the jet has explicit_ghosts (knowing that there is an +// area support) +bool BottomUpSoftDrop::_check_explicit_ghosts(const PseudoJet &jet) const { + // if the jet comes from a Clustering check explicit ghosts in that + // clustering + if (jet.has_associated_cluster_sequence()) + return jet.validated_csab()->has_explicit_ghosts(); + + // if the jet has pieces, recurse in the pieces + if (jet.has_pieces()){ + vector pieces = jet.pieces(); + for (unsigned int i=0;ijet_def().has_same_recombiner(jet_def_for_recombiner); + + // otherwise, assign it. + jet_def_for_recombiner = jet.validated_cs()->jet_def(); + assigned = true; + return true; + } + + // if the jet has pieces, recurse in the pieces + if (jet.has_pieces()){ + vector pieces = jet.pieces(); + if (pieces.size() == 0) return false; + for (unsigned int i=0;i kept(internal_hist.size(), true); + const vector &sd_rej = softdrop_recombiner.rejected(); + for (unsigned int i=0;i internal2input(internal_hist.size()); + for (unsigned int i=0; i