CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Combinatorics.cc
Go to the documentation of this file.
1 /*
2  * Combinatorics.cpp
3  *
4  * 03/04/2006 kasselmann@physik.rwth-aachen.de
5  * 19/08/2007 giffels@physik.rwth-aachen.de
6  *
7  */
8 //framework
10 
11 // Own header
13 
14 using namespace std;
15 
16 // **************************************************************************
17 // Class Constructor
18 // **************************************************************************
19 Combinatorics::Combinatorics(Int_t SetQuantity, Int_t SubsetQuantity) :
20 
21  m_SetQuantity(SetQuantity),
22  m_SubsetQuantity(SubsetQuantity)
23 {
24  // Get permutations
26 }
27 
28 
29 // **************************************************************************
30 // Class Destructor
31 // **************************************************************************
33 {}
34 
35 
36 // **************************************************************************
37 // Get subset permutations
38 // **************************************************************************
39 vector <vector <UInt_t> > Combinatorics::GetPermutations()
40 {
41  return m_Permutations;
42 }
43 
44 // **************************************************************************
45 // Calculate all subset permutations
46 // **************************************************************************
48 {
50  {
51  edm::LogWarning("Combinatorics") << "[Combinatorics] No valid choice of set or subset!" << endl;
52  return -1;
53  }
54 
55  Int_t* currentSubset = new Int_t[m_SubsetQuantity];
56  Int_t* currentMapping = new Int_t[m_SetQuantity];
57 
58  initial_subset(m_SubsetQuantity, currentSubset);
59  do
60  {
61  initial_permutation(m_SetQuantity, currentMapping);
62  do
63  {
64  for (UShort_t i = 0; i < m_SubsetQuantity; i++)
65  {
66  m_Subset.push_back(currentSubset[currentMapping[i]]);
67  }
68  m_Permutations.push_back(m_Subset);
69  m_Subset.clear();
70  }
71  while (next_permutation(m_SubsetQuantity, currentMapping));
72  }
73  while (next_subset(m_SetQuantity, m_SubsetQuantity, currentSubset));
74 
75  delete[] currentSubset;
76  delete[] currentMapping;
77 
78  return 0;
79 }
80 
81 
82 // **************************************************************************
83 // Build initial permutation
84 // **************************************************************************
85 void Combinatorics::initial_permutation(int size, int *permutation)
86 {
87  for (int i = 0; i < size; i++)
88  {
89  permutation[i] = i;
90  }
91 }
92 
93 
94 // **************************************************************************
95 // Build initial subset
96 // **************************************************************************
97 void Combinatorics::initial_subset(int k, int *subset)
98 {
99  for (int i = 0; i < k; i++)
100  {
101  subset[i] = i;
102  }
103 }
104 
105 
106 // **************************************************************************
107 // Get next permutation if a next permutation exists
108 // **************************************************************************
109 Bool_t Combinatorics::next_permutation(int size, int *permutation)
110 {
111  int i, j, k;
112  if (size < 2) return false;
113  i = size - 2;
114 
115  while ((permutation[i] > permutation[i+1]) && (i != 0))
116  {
117  i--;
118  }
119  if ((i == 0) && (permutation[0] > permutation[1])) return false;
120 
121  k = i + 1;
122  for (j = i + 2; j < size; j++ )
123  {
124  if ((permutation[j] > permutation[i]) && (permutation[j] < permutation[k]))
125  {
126  k = j;
127  }
128  }
129 
130  // swap i and k
131  {
132  int tmp = permutation[i];
133  permutation[i] = permutation[k];
134  permutation[k] = tmp;
135  }
136 
137  for (j = i + 1; j <= ((size + i) / 2); j++)
138  {
139  int tmp = permutation[j];
140  permutation[j] = permutation[size + i - j];
141  permutation[size + i - j] = tmp;
142  }
143 
144  // return whether a next permutation exists
145  return true;
146 }
147 
148 
149 // **************************************************************************
150 // Get next subset if a next subset exists
151 //
152 // n: the size of the set
153 // k: the size of the subset
154 // **************************************************************************
155 Bool_t Combinatorics::next_subset(int n, int k, int *subset)
156 {
157  int i;
158  int j;
159  int jsave;
160  bool done;
161 
162  if ( subset[0] < n-k )
163  {
164  done = false;
165  jsave = k-1;
166  for ( j = 0; j < k-1; j++ )
167  {
168  if ( subset[j] + 1 < subset[j+1] )
169  {
170  jsave = j;
171  break;
172  }
173  }
174  for ( i = 0; i < jsave; i++ ) subset[i] = i;
175  subset[jsave] = subset[jsave] + 1;
176  }
177  else
178  {
179  done = true;
180  }
181  // return whether a next subset exists
182  return !done;
183 }
184 
185 
186 // **************************************************************************
187 // Get all subset combinations
188 // **************************************************************************
189 vector < vector <UInt_t> > Combinatorics::GetCombinations()
190 {
191  if (m_Permutations.size() == 0)
192  {
193  LogDebug("Combinatorics") << "Nothing to do." << endl;
194  return m_Combinations;
195  }
196 
197  m_Combinations.push_back(m_Permutations.at(0));
198 
199  for (UInt_t i = 1; i < m_Permutations.size(); i++)
200  {
201  if (!EqualPermutation(m_Combinations.back(), m_Permutations.at(i)))
202  {
203  m_Combinations.push_back(m_Permutations.at(i));
204  }
205  }
206  return m_Combinations;
207 }
208 
209 
210 // **************************************************************************
211 // Returns true if two permutations of four "int" are equal
212 // (Equal means e.g.: 0123 = 1023 = 0132 = 1032)
213 // **************************************************************************
214 Int_t Combinatorics::EqualPermutation(const vector<UInt_t>& p1, const vector <UInt_t>& p2)
215 {
216  if (p1.size() != p2.size())
217  {
218  edm::LogWarning("Combinatorics") << "[Combinatorics::EqualPermutation] permutations have different size!" << endl;
219  return -1;
220  }
221 
222  Float_t p1_sum = 0.0;
223  Float_t p2_sum = 0.0;
224 
225  // Check whether permutations are equal (2^index)
226  for (UInt_t i = 0; i < p1.size(); i++) p1_sum += (1 << p1.at(i));
227  for (UInt_t i = 0; i < p2.size(); i++) p2_sum += (1 << p2.at(i));
228 
229  return (p1_sum == p2_sum ? 1 : 0);
230 }
231 
232 
233 // **************************************************************************
234 // Get combinations: 4 out of n
235 //
236 // (The order of the first and second two is not important!
237 // 0123 = 1023 = 0132 = 1032 are equal therefore)
238 // **************************************************************************
239 vector < vector <UInt_t> > Combinatorics::GetCombinations_2_2()
240 {
241  // combination vector returned
242  vector< vector <UInt_t> > FinalCombinations;
243 
244  if (m_Permutations.size() == 0)
245  {
246  LogDebug("Combinatorics") << "[Combinatorics::GetCombinations_2_2] Nothing to do." << endl;
247  return FinalCombinations;
248  }
249 
250  // So far only for subsets of four indices
251  if (m_SubsetQuantity != 4)
252  {
253  edm::LogWarning("Combinatorics") << "[Combinatorics::GetCombinations_2_2] Subset must be 4." << endl;
254  return FinalCombinations;
255  }
256 
257  // Skip specific permutations
258  Skip_2_2(m_Permutations, FinalCombinations);
259 
260  return FinalCombinations;
261 }
262 
263 
264 // **************************************************************************
265 // Get combinations: 4 out of n
266 //
267 // (The order of the last two is important only:
268 // 0123 = 1023 are equal therefore)
269 // **************************************************************************
270 vector < vector <UInt_t> > Combinatorics::GetCombinations_2_0()
271 {
272  // combination vector returned
273  vector< vector <UInt_t> > FinalCombinations;
274 
275  if (m_Permutations.size() == 0)
276  {
277  LogDebug("Combinatorics") << "[Combinatorics::GetCombinations_2_0] Nothing to do." << endl;
278  return FinalCombinations;
279  }
280 
281  // So far only for subsets of four indices
282  if (m_SubsetQuantity != 4)
283  {
284  edm::LogWarning("Combinatorics") << "[Combinatorics::GetCombinations_2_0] Subset must be 4." << endl;
285  return FinalCombinations;
286  }
287 
288  // Skip specific permutations
289  Skip_2_0(m_Permutations, FinalCombinations);
290 
291  return FinalCombinations;
292 }
293 
294 
295 // **************************************************************************
296 // Skip permutation from p1 if already existing in p2
297 // **************************************************************************
298 void Combinatorics::Skip_2_0(const vector <vector <UInt_t> >& p1,
299  vector <vector <UInt_t> > &p2)
300 {
301  Bool_t Skip = kFALSE;
302 
303  p2.push_back(p1.at(0));
304 
305  for (UShort_t i = 1; i < p1.size(); i++)
306  {
307  for (UShort_t j = 0; j < p2.size(); j++)
308  {
309  if (EqualPermutation_2_0(p1.at(i), p2.at(j)))
310  {
311  Skip = kTRUE;
312  }
313  }
314  if (!Skip) p2.push_back(p1.at(i));
315 
316  Skip = kFALSE;
317  }
318 }
319 
320 
321 // **************************************************************************
322 // Skip permutation from p1 if already existing in p2
323 // **************************************************************************
324 void Combinatorics::Skip_2_2(const vector <vector <UInt_t> >& p1,
325  vector <vector <UInt_t> > &p2)
326 {
327  Bool_t Skip = kFALSE;
328 
329  p2.push_back(p1.at(0));
330 
331  for (UShort_t i = 1; i < p1.size(); i++)
332  {
333  for (UShort_t j = 0; j < p2.size(); j++)
334  {
335  if (EqualPermutation_2_2(p1.at(i), p2.at(j)))
336  {
337  Skip = kTRUE;
338  }
339  }
340  if (!Skip) p2.push_back(p1.at(i));
341 
342  Skip = kFALSE;
343  }
344 }
345 
346 
347 // **************************************************************************
348 // Returns true if the two first digm_ of two permutations are equal
349 // e.g.: 0123 = 1023
350 // **************************************************************************
351 Int_t Combinatorics::EqualPermutation_2_0(const vector<UInt_t>& p1, const vector<UInt_t>& p2)
352 {
353  if (p1.size() < 2)
354  {
355  edm::LogWarning("Combinatorics") << "[Combinatorics::EqualPermutation_2_0] permutation has wrong size!" << endl;
356  return -1;
357  }
358 
359  // Check whether permutations are equal
360  if ( ((1 << p1.at(0)) + (1 << p1.at(1)) == (1 << p2.at(0)) + (1 << p2.at(1)) ) &&
361  p1.at(2) == p2.at(2) && p1.at(3) == p2.at(3) )
362  {
363  return 1;
364  }
365  return 0;
366 }
367 
368 
369 // **************************************************************************
370 // Returns true if two permutations of four "int" are equal
371 // e.g.: 0123 = 1023 = 0132 = 1032
372 // **************************************************************************
373 Int_t Combinatorics::EqualPermutation_2_2(const vector<UInt_t>& p1, const vector<UInt_t>& p2)
374 {
375  // Returns true if two permutations of four "int" are equal
376  // (equal means e.g.: 0123 = 1023 = 0132 = 1032)
377 
378  if (p1.size() != 4 && p2.size() != 4)
379  {
380  edm::LogWarning("Combinatorics") << "[Combinatorics::EqualPermutationTwoByTwo] permutation(s) have wrong size!" << endl;
381  return -1;
382  }
383 
384  // Check whether permutations are equal (2^index)
385  if ( ((1 << p1.at(0)) + (1 << p1.at(1)) == (1 << p2.at(0)) + (1 << p2.at(1)) ) &&
386  ((1 << p1.at(2)) + (1 << p1.at(3)) == (1 << p2.at(2)) + (1 << p2.at(3)) ) )
387  {
388  return 1;
389  }
390  return 0;
391 }
392 
393 
394 // **************************************************************************
395 // Returns true if two permutations of four are "equal"
396 // e.g.: 0123 = 1023
397 // **************************************************************************
398 Int_t Combinatorics::EqualPermutation_N_1(const vector<UInt_t>& p1,const vector<UInt_t>& p2)
399 {
400  // Returns true if two permutations of four "int" are equal
401  // (equal means e.g.: 012 = 102)
402 
403  if (p1.size() != p2.size())
404  {
405  edm::LogWarning("Combinatorics") << "[Combinatorics::EqualPermutationTwoByTwo] permutation(s) have wrong size!" << endl;
406  return -1;
407  }
408 
409  return (EqualPermutation(p1, p2) && p1.back() == p2.back() ? 1 : 0);
410 }
411 
412 
413 // **************************************************************************
414 // Get combinations "N by 1"
415 // **************************************************************************
416 vector < vector <UInt_t> > Combinatorics::GetCombinations_N_1()
417 {
418  // Get combinations
419  m_Combinations.clear();
420  GetCombinations();
421 
422  // combination vector returned
423  vector < vector <UInt_t> > FinalCombinations;
424 
425  if (m_Combinations.size() == 0)
426  {
427  LogDebug("Combinatorics") << "[Combinatorics::GetCombinationsThreeByOne] Nothing to do." << endl;
428  return FinalCombinations;
429  }
430 
431  for (UInt_t i = 0; i < m_Combinations.size(); i++)
432  {
433  vector <UInt_t> RotatingPermutation = m_Combinations.at(i);
434  FinalCombinations.push_back(m_Combinations.at(i));
435 
436  for (UInt_t j = 1; j < RotatingPermutation.size(); j++)
437  {
438  FinalCombinations.push_back(Rotate(RotatingPermutation,j));
439  }
440  }
441  return FinalCombinations;
442 }
443 
444 
445 // **************************************************************************
446 // Rotate permutation to the "left" by n digm_
447 // **************************************************************************
448 vector<UInt_t> Combinatorics::Rotate(const vector <UInt_t>& permutation, UInt_t digm_)
449 {
450  vector<UInt_t> p;
451  vector<UInt_t> tmp;
452 
453  if (permutation.size() <= digm_)
454  {
455  edm::LogWarning("Combinatorics") << "[Combinatorics::Rotate] WARNING: More rotations than digm_ in permutation!" << endl;
456  }
457 
458  // Save the first i digm_
459  for (UInt_t i = 0; i < digm_; i++)
460  {
461  tmp.push_back(permutation.at(i));
462  }
463  for (UInt_t j = 0; j < permutation.size() - digm_; j++)
464  {
465  p.push_back(permutation.at(j + digm_));
466  }
467  for (UInt_t k = 0; k < digm_; k++) p.push_back(tmp.at(k));
468 
469  return p;
470 }
471 
472 
473 // **************************************************************************
474 // Print one permutation
475 // **************************************************************************
476 void Combinatorics::Print(const vector<UInt_t>& p)
477 {
478  // Print permutations
479  for (UShort_t i = 0; i < p.size(); i++)
480  {
481  LogDebug("Combinatorics") << (p.at(i));
482  }
483  LogDebug("Combinatorics") << endl;
484 }
485 
486 
487 // **************************************************************************
488 // Print permutations
489 // **************************************************************************
490 void Combinatorics::Print(const vector <vector <UInt_t> >& p)
491 {
492  LogDebug("Combinatorics") << "**************" << endl;
493  LogDebug("Combinatorics") << "Permutations: " << p.size() << endl;
494 
495  // Print permutations
496  for (UShort_t i = 0; i < p.size(); i++)
497  {
498  for(UShort_t j = 0; j < (p.at(0)).size(); j++) LogDebug("Combinatorics") << (p.at(i)).at(j);
499  LogDebug("Combinatorics") << endl;
500  }
501 }
502 
503 
504 
505 
506 
507 
508 
509 
510 
511 
512 
#define LogDebug(id)
int i
Definition: DBlmapReader.cc:9
Int_t EqualPermutation_N_1(const std::vector< UInt_t > &permutation1, const std::vector< UInt_t > &permutation2)
std::vector< std::vector< UInt_t > > GetCombinations()
Int_t EqualPermutation_2_2(const std::vector< UInt_t > &permutation1, const std::vector< UInt_t > &permutation2)
std::vector< std::vector< UInt_t > > m_Permutations
Definition: Combinatorics.h:63
std::vector< std::vector< UInt_t > > GetCombinations_N_1()
Combinatorics(Int_t Set, Int_t Subset)
Int_t CalculatePermutations()
const Int_t m_SubsetQuantity
Definition: Combinatorics.h:60
std::vector< UInt_t > m_Subset
Definition: Combinatorics.h:62
int j
Definition: DBlmapReader.cc:9
Bool_t next_subset(int n, int k, int *subset)
const Int_t m_SetQuantity
Definition: Combinatorics.h:59
std::vector< std::vector< UInt_t > > GetCombinations_2_0()
Int_t EqualPermutation_2_0(const std::vector< UInt_t > &permutation1, const std::vector< UInt_t > &permutation2)
double p2[4]
Definition: TauolaWrapper.h:90
virtual ~Combinatorics()
std::vector< std::vector< UInt_t > > m_Combinations
Definition: Combinatorics.h:64
Bool_t next_permutation(int size, int *permutation)
void Skip_2_0(const std::vector< std::vector< UInt_t > > &permutation1, std::vector< std::vector< UInt_t > > &permutation2)
std::vector< std::vector< UInt_t > > GetPermutations()
std::vector< std::vector< UInt_t > > GetCombinations_2_2()
void Print(const std::vector< UInt_t > &permutation)
Int_t EqualPermutation(const std::vector< UInt_t > &permutation1, const std::vector< UInt_t > &permutation2)
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
void initial_permutation(int size, int *permutation)
double p1[4]
Definition: TauolaWrapper.h:89
void initial_subset(int k, int *subset)
tuple size
Write out results.
std::vector< UInt_t > Rotate(const std::vector< UInt_t > &permutation, UInt_t digits)
void Skip_2_2(const std::vector< std::vector< UInt_t > > &permutation1, std::vector< std::vector< UInt_t > > &permutation2)