Fork me on GitHub

source: git/external/fastjet/plugins/SISCone/hash.cc@ 7bb13cd

ImprovedOutputFile Timing dual_readout llp
Last change on this file since 7bb13cd was 35cdc46, checked in by Pavel Demin <demin@…>, 10 years ago

upgrade FastJet to version 3.1.0-beta.1, upgrade Nsubjettiness to version 2.1.0, add SoftKiller version 1.0.0

  • Property mode set to 100644
File size: 7.5 KB
Line 
1///////////////////////////////////////////////////////////////////////////////
2// File: hash.cpp //
3// Description: source file for classes hash_element and hash_cones //
4// This file is part of the SISCone project. //
5// For more details, see http://projects.hepforge.org/siscone //
6// //
7// Copyright (c) 2006 Gavin Salam and Gregory Soyez //
8// //
9// This program is free software; you can redistribute it and/or modify //
10// it under the terms of the GNU General Public License as published by //
11// the Free Software Foundation; either version 2 of the License, or //
12// (at your option) any later version. //
13// //
14// This program is distributed in the hope that it will be useful, //
15// but WITHOUT ANY WARRANTY; without even the implied warranty of //
16// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the //
17// GNU General Public License for more details. //
18// //
19// You should have received a copy of the GNU General Public License //
20// along with this program; if not, write to the Free Software //
21// Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA //
22// //
23// $Revision:: 225 $//
24// $Date:: 2008-05-20 16:59:47 +0200 (Tue, 20 May 2008) $//
25///////////////////////////////////////////////////////////////////////////////
26
27#include <math.h>
28#include <stdio.h>
29#include "hash.h"
30#include <iostream>
31
32namespace siscone{
33
34using namespace std;
35
36/**************************************************************
37 * implementation of hash_cones *
38 * list of cones candidates. *
39 * We store in this class all the hash_elements and give *
40 * functions to manipulate them. *
41 **************************************************************/
42
43// constructor with initialisation
44// - _Np number of particles
45// - _R2 cone radius (squared)
46//-----------------------------------
47hash_cones::hash_cones(int _Np, double _R2){
48 int i;
49
50 n_cones = 0;
51#ifdef DEBUG_STABLE_CONES
52 n_occupied_cells = 0;
53#endif
54
55 // determine hash size
56 // for a ymax=5 and R=0.7, we observed an occupancy around 1/8 N^2 ~ N2 R2/4
57 //mask = 1 << (int) (2*log(double(_Np))/log(2.0));
58 //if (mask<=1) mask=2;
59 int nbits = (int) (log(_Np*_R2*_Np/4.0)/log(2.0));
60 if (nbits<1) nbits=1;
61 mask = 1 << nbits;
62
63 // create hash
64 hash_array = new hash_element*[mask];
65 mask--;
66
67 // set the array to 0
68 //? needed ?
69 for (i=0;i<mask+1;i++)
70 hash_array[i] = NULL;
71
72 R2 = _R2;
73}
74
75// destructor
76//------------
77hash_cones::~hash_cones(){
78 int i;
79 hash_element *elm;
80
81 for (i=0;i<mask+1;i++){
82 while (hash_array[i]!=NULL){
83 elm = hash_array[i];
84 hash_array[i] = hash_array[i]->next;
85 delete elm;
86 }
87 }
88
89 delete[] hash_array;
90}
91
92
93/*
94 * insert a new candidate into the hash.
95 * - v 4-momentum of the cone to add
96 * - parent parent particle defining the cone
97 * - child child particle defining the cone
98 * - p_io whether the parent has to belong to the cone or not
99 * - c_io whether the child has to belong to the cone or not
100 * return 0 on success, 1 on error
101 ***********************************************************************/
102int hash_cones::insert(Cmomentum *v, Cmomentum *parent, Cmomentum *child, bool p_io, bool c_io){
103 hash_element *elm;
104 int index = (v->ref.ref[0]) & mask;
105
106 // check the array cell corresponding to our reference
107 elm = hash_array[index];
108
109#ifdef DEBUG_STABLE_CONES
110 if (elm==NULL)
111 n_occupied_cells++;
112#endif
113
114 do{
115 // if it is not present, add it
116 if (elm==NULL){
117 // create element
118 elm = new hash_element;
119
120 // set its varibles
121 // Note: at this level, eta and phi have already been computed
122 // through Cmomentum::build_etaphi.
123 elm->ref = v->ref;
124
125 //compute vectors centre
126 v->build_etaphi();
127 elm->eta = v->eta;
128 elm->phi = v->phi;
129 // if at least one of the two is_inside tests gives a result != from the expected,
130 // the || will be true hence !(...) false as wanted
131 elm->is_stable = !((is_inside(v, parent)^p_io)||(is_inside(v, child)^c_io));
132 //cout << "-- new status of " << v->ref[0] << ":" << elm->is_stable << endl;
133
134 // update hash
135 elm->next = hash_array[index];
136 hash_array[index] = elm;
137
138 n_cones++;
139 return 0;
140 }
141
142 // if the cone is already there, simply update stability status
143 if (v->ref == elm->ref){
144 // there is only an update to perform to see if the cone is still stable
145 if (elm->is_stable){
146 v->build_etaphi();
147 elm->is_stable = !((is_inside(v, parent)^p_io)||(is_inside(v, child)^c_io));
148 //cout << " parent/child: "
149 // << parent->ref[0] << ":" << is_inside(v, parent) << ":" << p_io << " "
150 // << child->ref[0] << ":" << is_inside(v, child) << ":" << c_io << endl;
151 //cout << "-- rep status of " << v->ref[0] << ":" << elm->is_stable << endl;
152 //cout << v->eta << " " << v->phi << endl;
153 //cout << (child->eta) << " " << child->phi << endl;
154 }
155 return 0;
156 }
157
158 elm = elm->next;
159 } while (1);
160
161 return 1;
162}
163
164/*
165 * insert a new candidate into the hash.
166 * - v 4-momentum of te cone to add
167 * Note, in this case, we assume stability. We also assume
168 * that eta and phi are computed for v
169 * return 0 on success, 1 on error
170 ***********************************************************************/
171int hash_cones::insert(Cmomentum *v){
172 hash_element *elm;
173 int index = (v->ref.ref[0]) & mask;
174 //cout << "-- stable candidate: " << v->ref[0] << ":" << endl;
175
176 // check the array cell corresponding to our reference
177 elm = hash_array[index];
178 do{
179 // if it is not present, add it
180 if (elm==NULL){
181 // create element
182 elm = new hash_element;
183
184 // set its varibles
185 // Note: at this level, eta and phi have already been computed
186 // through Cmomentum::build_etaphi.
187 elm->ref = v->ref;
188 elm->eta = v->eta;
189 elm->phi = v->phi;
190 elm->is_stable = true;
191
192 // update hash
193 elm->next = hash_array[index];
194 hash_array[index] = elm;
195
196 n_cones++;
197 return 0;
198 }
199
200 // if the cone is already there, we have nothing to do
201 if (v->ref == elm->ref){
202 return 0;
203 }
204
205 elm = elm->next;
206 } while (1);
207
208 return 1;
209}
210
211/*
212 * test if a particle is inside a cone of given centre.
213 * check if the particle of coordinates 'v' is inside the circle of radius R
214 * centered at 'centre'.
215 * - centre centre of the circle
216 * - v particle to test
217 * return true if inside, false if outside
218 ******************************************************************************/
219inline bool hash_cones::is_inside(Cmomentum *centre, Cmomentum *v){
220 double dx, dy;
221
222 dx = centre->eta - v->eta;
223 dy = fabs(centre->phi - v->phi);
224 if (dy>M_PI)
225 dy -= 2.0*M_PI;
226
227 return dx*dx+dy*dy<R2;
228}
229
230}
Note: See TracBrowser for help on using the repository browser.