Fork me on GitHub

source: git/external/fastjet/contribs/Nsubjettiness/Njettiness.cc@ 9e77ec8

ImprovedOutputFile Timing dual_readout llp
Last change on this file since 9e77ec8 was 1d208a2, checked in by Pavel Demin <pavel.demin@…>, 8 years ago

update FastJet library to 3.2.1 and Nsubjettiness library to 2.2.4

  • Property mode set to 100644
File size: 8.3 KB
RevLine 
[9687203]1// Nsubjettiness Package
2// Questions/Comments? jthaler@jthaler.net
3//
4// Copyright (c) 2011-14
5// Jesse Thaler, Ken Van Tilburg, Christopher K. Vermilion, and TJ Wilkason
6//
[1d208a2]7// $Id: Njettiness.cc 933 2016-04-04 22:23:32Z jthaler $
[9687203]8//----------------------------------------------------------------------
9// This file is part of FastJet contrib.
10//
11// It is free software; you can redistribute it and/or modify it under
12// the terms of the GNU General Public License as published by the
13// Free Software Foundation; either version 2 of the License, or (at
14// your option) any later version.
15//
16// It is distributed in the hope that it will be useful, but WITHOUT
17// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
18// or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
19// License for more details.
20//
21// You should have received a copy of the GNU General Public License
22// along with this code. If not, see <http://www.gnu.org/licenses/>.
23//----------------------------------------------------------------------
24
25#include "Njettiness.hh"
26
27FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
28
29namespace contrib {
30
31
32///////
33//
34// Main Njettiness Class
35//
36///////
37
[973b92a]38LimitedWarning Njettiness::_old_measure_warning;
39LimitedWarning Njettiness::_old_axes_warning;
40
41
42// Constructor
[35cdc46]43Njettiness::Njettiness(const AxesDefinition & axes_def, const MeasureDefinition & measure_def)
[973b92a]44: _axes_def(axes_def.create()), _measure_def(measure_def.create()) {}
45
46// setAxes for Manual mode
47void Njettiness::setAxes(const std::vector<fastjet::PseudoJet> & myAxes) {
[1d208a2]48 if (_axes_def->needsManualAxes()) {
[973b92a]49 _currentAxes = myAxes;
50 } else {
51 throw Error("You can only use setAxes for manual AxesDefinitions");
52 }
[9687203]53}
[973b92a]54
55// Calculates and returns all TauComponents that user would want.
56// This information is stored in _current_tau_components for later access as well.
57TauComponents Njettiness::getTauComponents(unsigned n_jets, const std::vector<fastjet::PseudoJet> & inputJets) const {
58 if (inputJets.size() <= n_jets) { //if not enough particles, return zero
59 _currentAxes = inputJets;
60 _currentAxes.resize(n_jets,fastjet::PseudoJet(0.0,0.0,0.0,0.0));
[1d208a2]61
62 // Put in empty tau components
63 std::vector<double> dummy_jet_pieces;
64 _current_tau_components = TauComponents(UNDEFINED_SHAPE,
65 dummy_jet_pieces,
66 0.0,
67 1.0,
68 _currentAxes,
69 _currentAxes
70 );
[973b92a]71 _seedAxes = _currentAxes;
72 _currentPartition = TauPartition(n_jets); // empty partition
73 } else {
[1d208a2]74 assert(_axes_def); // this should never fail.
[973b92a]75
[1d208a2]76 if (_axes_def->needsManualAxes()) { // if manual mode
[973b92a]77 // take current axes as seeds
78 _seedAxes = _currentAxes;
79
80 // refine axes if requested
[1d208a2]81 _currentAxes = _axes_def->get_refined_axes(n_jets,inputJets,_seedAxes, _measure_def.get());
[973b92a]82 } else { // non-manual axes
83
84 //set starting point for minimization
[1d208a2]85 _seedAxes = _axes_def->get_starting_axes(n_jets,inputJets,_measure_def.get());
[973b92a]86
87 // refine axes as needed
[1d208a2]88 _currentAxes = _axes_def->get_refined_axes(n_jets,inputJets,_seedAxes, _measure_def.get());
[973b92a]89
90 // NOTE: The above two function calls are combined in "AxesDefinition::get_axes"
91 // but are separated here to allow seed axes to be stored.
92 }
93
94 // Find and store partition
95 _currentPartition = _measure_def->get_partition(inputJets,_currentAxes);
96
97 // Find and store tau value
98 _current_tau_components = _measure_def->component_result_from_partition(_currentPartition, _currentAxes); // sets current Tau Values
99 }
100 return _current_tau_components;
[35cdc46]101}
102
[973b92a]103
104///////
105//
106// Below is code for backward compatibility to use the old interface.
107// May be deleted in a future version
108//
109///////
110
111Njettiness::Njettiness(AxesMode axes_mode, const MeasureDefinition & measure_def)
112: _axes_def(createAxesDef(axes_mode)), _measure_def(measure_def.create()) {}
113
[35cdc46]114// Convert from MeasureMode enum to MeasureDefinition
115// This returns a pointer that will be claimed by a SharedPtr
116MeasureDefinition* Njettiness::createMeasureDef(MeasureMode measure_mode, int num_para, double para1, double para2, double para3) const {
[973b92a]117
118 _old_measure_warning.warn("Njettiness::createMeasureDef: You are using the old MeasureMode way of specifying N-subjettiness measures. This is deprecated as of v2.1 and will be removed in v3.0. Please use MeasureDefinition instead.");
119
[9687203]120 // definition of maximum Rcutoff for non-cutoff measures, changed later by other measures
121 double Rcutoff = std::numeric_limits<double>::max(); //large number
[35cdc46]122 // Most (but not all) measures have some kind of beta value
123 double beta = std::numeric_limits<double>::quiet_NaN();
[9687203]124 // The normalized measures have an R0 value.
[35cdc46]125 double R0 = std::numeric_limits<double>::quiet_NaN();
126
[9687203]127 // Find the MeasureFunction and set the parameters.
128 switch (measure_mode) {
129 case normalized_measure:
130 beta = para1;
131 R0 = para2;
[35cdc46]132 if(num_para == 2) {
133 return new NormalizedMeasure(beta,R0);
134 } else {
135 throw Error("normalized_measure needs 2 parameters (beta and R0)");
136 }
[9687203]137 break;
138 case unnormalized_measure:
139 beta = para1;
[35cdc46]140 if(num_para == 1) {
141 return new UnnormalizedMeasure(beta);
142 } else {
143 throw Error("unnormalized_measure needs 1 parameter (beta)");
144 }
[9687203]145 break;
146 case geometric_measure:
[973b92a]147 throw Error("This class has been removed. Please use OriginalGeometricMeasure, ModifiedGeometricMeasure, or ConicalGeometricMeasure with the new Njettiness constructor.");
[9687203]148 break;
149 case normalized_cutoff_measure:
150 beta = para1;
151 R0 = para2;
152 Rcutoff = para3; //Rcutoff parameter is 3rd parameter in normalized_cutoff_measure
[35cdc46]153 if (num_para == 3) {
154 return new NormalizedCutoffMeasure(beta,R0,Rcutoff);
155 } else {
156 throw Error("normalized_cutoff_measure has 3 parameters (beta, R0, Rcutoff)");
157 }
[9687203]158 break;
159 case unnormalized_cutoff_measure:
160 beta = para1;
161 Rcutoff = para2; //Rcutoff parameter is 2nd parameter in normalized_cutoff_measure
[35cdc46]162 if (num_para == 2) {
163 return new UnnormalizedCutoffMeasure(beta,Rcutoff);
164 } else {
165 throw Error("unnormalized_cutoff_measure has 2 parameters (beta, Rcutoff)");
166 }
[9687203]167 break;
168 case geometric_cutoff_measure:
[973b92a]169 throw Error("This class has been removed. Please use OriginalGeometricMeasure, ModifiedGeometricMeasure, or ConicalGeometricMeasure with the new Njettiness constructor.");
[9687203]170 default:
171 assert(false);
172 break;
[35cdc46]173 }
174 return NULL;
175}
[9687203]176
[35cdc46]177// Convert from AxesMode enum to AxesDefinition
178// This returns a pointer that will be claimed by a SharedPtr
179AxesDefinition* Njettiness::createAxesDef(Njettiness::AxesMode axes_mode) const {
180
[973b92a]181 _old_axes_warning.warn("Njettiness::createAxesDef: You are using the old AxesMode way of specifying N-subjettiness axes. This is deprecated as of v2.1 and will be removed in v3.0. Please use AxesDefinition instead.");
182
183
[9687203]184 switch (axes_mode) {
185 case wta_kt_axes:
[35cdc46]186 return new WTA_KT_Axes();
[9687203]187 case wta_ca_axes:
[35cdc46]188 return new WTA_CA_Axes();
[9687203]189 case kt_axes:
[35cdc46]190 return new KT_Axes();
[9687203]191 case ca_axes:
[35cdc46]192 return new CA_Axes();
[9687203]193 case antikt_0p2_axes:
[35cdc46]194 return new AntiKT_Axes(0.2);
[9687203]195 case onepass_wta_kt_axes:
[35cdc46]196 return new OnePass_WTA_KT_Axes();
[9687203]197 case onepass_wta_ca_axes:
[35cdc46]198 return new OnePass_WTA_CA_Axes();
[9687203]199 case onepass_kt_axes:
[35cdc46]200 return new OnePass_KT_Axes();
[9687203]201 case onepass_ca_axes:
[35cdc46]202 return new OnePass_CA_Axes();
[9687203]203 case onepass_antikt_0p2_axes:
[35cdc46]204 return new OnePass_AntiKT_Axes(0.2);
[9687203]205 case onepass_manual_axes:
[35cdc46]206 return new OnePass_Manual_Axes();
207 case min_axes:
208 return new MultiPass_Axes(100);
[9687203]209 case manual_axes:
[35cdc46]210 return new Manual_Axes();
[9687203]211 default:
212 assert(false);
[35cdc46]213 return NULL;
214 }
215}
[9687203]216
[35cdc46]217
[9687203]218
219
220
221} // namespace contrib
222
223FASTJET_END_NAMESPACE
Note: See TracBrowser for help on using the repository browser.