ROOT logo

From $ROOTSYS/tutorials/roostats/rs_numberCountingCombination.C

/////////////////////////////////////////////////////////////////////////
//
// 'Number Counting Example' RooStats tutorial macro #100
// author: Kyle Cranmer
// date Nov. 2008 
//
// This tutorial shows an example of a combination of 
// two searches using number counting with background uncertainty.
//
// The macro uses a RooStats "factory" to construct a PDF
// that represents the two number counting analyses with background 
// uncertainties.  The uncertainties are taken into account by 
// considering a sideband measurement of a size that corresponds to the
// background uncertainty.  The problem has been studied in these references:
//   http://arxiv.org/abs/physics/0511028
//   http://arxiv.org/abs/physics/0702156
//   http://cdsweb.cern.ch/record/1099969?ln=en
//
// After using the factory to make the model, we use a RooStats 
// ProfileLikelihoodCalculator for a Hypothesis test and a confidence interval.
// The calculator takes into account systematics by eliminating nuisance parameters
// with the profile likelihood.  This is equivalent to the method of MINOS.
//
/////////////////////////////////////////////////////////////////////////

#ifndef __CINT__
#include "RooGlobalFunc.h"
#endif
#include "RooStats/ProfileLikelihoodCalculator.h"
#include "RooStats/NumberCountingPdfFactory.h"
#include "RooStats/ConfInterval.h"
#include "RooStats/HypoTestResult.h"
#include "RooStats/LikelihoodIntervalPlot.h"
#include "RooRealVar.h"

// use this order for safety on library loading
using namespace RooFit ;
using namespace RooStats ;


// declare three variations on the same tutorial
void rs_numberCountingCombination_expected();
void rs_numberCountingCombination_observed();
void rs_numberCountingCombination_observedWithTau();

////////////////////////////////////////////
// main driver to choose one
void rs_numberCountingCombination(int flag=1)
{
  if(flag==1) 
    rs_numberCountingCombination_expected();
  if(flag==2) 
    rs_numberCountingCombination_observed();
  if(flag==3) 
    rs_numberCountingCombination_observedWithTau();
}

/////////////////////////////////////////////
void rs_numberCountingCombination_expected()
{

  /////////////////////////////////////////
  // An example of a number counting combination with two channels.
  // We consider both hypothesis testing and the equivalent confidence interval.
  /////////////////////////////////////////


  /////////////////////////////////////////
  // The Model building stage
  /////////////////////////////////////////

  // Step 1, define arrays with signal & bkg expectations and background uncertainties
  Double_t s[2] = {20.,10.};           // expected signal
  Double_t b[2] = {100.,100.};         // expected background
  Double_t db[2] = {.0100,.0100};      // fractional background uncertainty

  
  // Step 2, use a RooStats factory to build a PDF for a 
  // number counting combination and add it to the workspace.
  // We need to give the signal expectation to relate the masterSignal
  // to the signal contribution in the individual channels.
  // The model neglects correlations in background uncertainty, 
  // but they could be added without much change to the example.
  NumberCountingPdfFactory f;
  RooWorkspace* wspace = new RooWorkspace();
  f.AddModel(s,2,wspace,"TopLevelPdf", "masterSignal"); 

  // Step 3, use a RooStats factory to add datasets to the workspace.
  // Step 3a.
  // Add the expected data to the workspace
  f.AddExpData(s, b, db, 2, wspace, "ExpectedNumberCountingData");

  // see below for a printout of the workspace
  //  wspace->Print();  //uncomment to see structure of workspace

  /////////////////////////////////////////
  // The Hypothesis testing stage:
  /////////////////////////////////////////
  // Step 4, Define the null hypothesis for the calculator
  // Here you need to know the name of the variables corresponding to hypothesis.
  RooRealVar* mu = wspace->var("masterSignal"); 
  RooArgSet* poi = new RooArgSet(*mu); 
  RooArgSet* nullParams = new RooArgSet("nullParams");
  nullParams->addClone(*mu);
  // here we explicitly set the value of the parameters for the null
  nullParams->setRealValue("masterSignal",0); 

  // Step 5, Create a calculator for doing the hypothesis test.
  // because this is a 
  ProfileLikelihoodCalculator plc( *wspace->data("ExpectedNumberCountingData"),
				   *wspace->pdf("TopLevelPdf"), *poi, 0.05, nullParams);
				  

  // Step 6, Use the Calculator to get a HypoTestResult
  HypoTestResult* htr = plc.GetHypoTest();
  assert(htr != 0);
  cout << "-------------------------------------------------" << endl;
  cout << "The p-value for the null is " << htr->NullPValue() << endl;
  cout << "Corresponding to a signifcance of " << htr->Significance() << endl;
  cout << "-------------------------------------------------\n\n" << endl;

  /* expected case should return:
     -------------------------------------------------
     The p-value for the null is 0.015294
     Corresponding to a signifcance of 2.16239
     -------------------------------------------------
  */

  //////////////////////////////////////////
  // Confidence Interval Stage

  // Step 8, Here we re-use the ProfileLikelihoodCalculator to return a confidence interval.
  // We need to specify what are our parameters of interest
  RooArgSet* paramsOfInterest = nullParams; // they are the same as before in this case
  plc.SetParameters(*paramsOfInterest);
  LikelihoodInterval* lrint = (LikelihoodInterval*) plc.GetInterval();  // that was easy.
  lrint->SetConfidenceLevel(0.95);

  // Step 9, make a plot of the likelihood ratio and the interval obtained
  //paramsOfInterest->setRealValue("masterSignal",1.); 
  // find limits
  double lower = lrint->LowerLimit(*mu);
  double upper = lrint->UpperLimit(*mu);

  LikelihoodIntervalPlot lrPlot(lrint);
  lrPlot.SetMaximum(3.);
  lrPlot.Draw();

  // Step 10a. Get upper and lower limits
  cout << "lower limit on master signal = " <<  lower << endl;
  cout << "upper limit on master signal = " <<  upper << endl;


  // Step 10b, Ask if masterSignal=0 is in the interval.
  // Note, this is equivalent to the question of a 2-sigma hypothesis test: 
  // "is the parameter point masterSignal=0 inside the 95% confidence interval?"
  // Since the signficance of the Hypothesis test was > 2-sigma it should not be: 
  // eg. we exclude masterSignal=0 at 95% confidence.
  paramsOfInterest->setRealValue("masterSignal",0.); 
  cout << "-------------------------------------------------" << endl;
  std::cout << "Consider this parameter point:" << std::endl;
  paramsOfInterest->first()->Print();
  if( lrint->IsInInterval(*paramsOfInterest) )
    std::cout << "It IS in the interval."  << std::endl;
  else
    std::cout << "It is NOT in the interval."  << std::endl;
  cout << "-------------------------------------------------\n\n" << endl;

  // Step 10c, We also ask about the parameter point masterSignal=2, which is inside the interval.
  paramsOfInterest->setRealValue("masterSignal",2.); 
  cout << "-------------------------------------------------" << endl;
  std::cout << "Consider this parameter point:" << std::endl;
  paramsOfInterest->first()->Print();
  if( lrint->IsInInterval(*paramsOfInterest) )
    std::cout << "It IS in the interval."  << std::endl;
  else
    std::cout << "It is NOT in the interval."  << std::endl;
  cout << "-------------------------------------------------\n\n" << endl;
  

  delete lrint;
  delete htr;
  delete wspace;
  delete poi; 
  delete nullParams;



  /*
  // Here's an example of what is in the workspace 
  //  wspace->Print();
  RooWorkspace(NumberCountingWS) Number Counting WS contents

  variables
  ---------
  (x_0,masterSignal,expected_s_0,b_0,y_0,tau_0,x_1,expected_s_1,b_1,y_1,tau_1)
  
  p.d.f.s
  -------
  RooProdPdf::joint[ pdfs=(sigRegion_0,sideband_0,sigRegion_1,sideband_1) ] = 2.20148e-08
  RooPoisson::sigRegion_0[ x=x_0 mean=splusb_0 ] = 0.036393
  RooPoisson::sideband_0[ x=y_0 mean=bTau_0 ] = 0.00398939
  RooPoisson::sigRegion_1[ x=x_1 mean=splusb_1 ] = 0.0380088
  RooPoisson::sideband_1[ x=y_1 mean=bTau_1 ] = 0.00398939
  
  functions
  --------
  RooAddition::splusb_0[ set1=(s_0,b_0) set2=() ] = 120
  RooProduct::s_0[ compRSet=(masterSignal,expected_s_0) compCSet=() ] = 20
  RooProduct::bTau_0[ compRSet=(b_0,tau_0) compCSet=() ] = 10000
  RooAddition::splusb_1[ set1=(s_1,b_1) set2=() ] = 110
  RooProduct::s_1[ compRSet=(masterSignal,expected_s_1) compCSet=() ] = 10
  RooProduct::bTau_1[ compRSet=(b_1,tau_1) compCSet=() ] = 10000
  
  datasets
  --------
  RooDataSet::ExpectedNumberCountingData(x_0,y_0,x_1,y_1)
  
  embedded precalculated expensive components
  -------------------------------------------
  */
  
}



void rs_numberCountingCombination_observed()
{

  /////////////////////////////////////////
  // The same example with observed data in a main
  // measurement and an background-only auxiliary 
  // measurement with a factor tau more background
  // than in the main measurement.

  /////////////////////////////////////////
  // The Model building stage
  /////////////////////////////////////////

  // Step 1, define arrays with signal & bkg expectations and background uncertainties
  // We still need the expectation to relate signal in different channels with the master signal
  Double_t s[2] = {20.,10.};           // expected signal

  
  // Step 2, use a RooStats factory to build a PDF for a 
  // number counting combination and add it to the workspace.
  // We need to give the signal expectation to relate the masterSignal
  // to the signal contribution in the individual channels.
  // The model neglects correlations in background uncertainty, 
  // but they could be added without much change to the example.
  NumberCountingPdfFactory f;
  RooWorkspace* wspace = new RooWorkspace();
  f.AddModel(s,2,wspace,"TopLevelPdf", "masterSignal"); 

  // Step 3, use a RooStats factory to add datasets to the workspace.
  // Add the observed data to the workspace
  Double_t mainMeas[2] = {123.,117.};      // observed main measurement
  Double_t bkgMeas[2] = {111.23,98.76};    // observed background
  Double_t dbMeas[2] = {.011,.0095};       // observed fractional background uncertainty
  f.AddData(mainMeas, bkgMeas, dbMeas, 2, wspace,"ObservedNumberCountingData");

  // see below for a printout of the workspace
  //  wspace->Print();  //uncomment to see structure of workspace

  /////////////////////////////////////////
  // The Hypothesis testing stage:
  /////////////////////////////////////////
  // Step 4, Define the null hypothesis for the calculator
  // Here you need to know the name of the variables corresponding to hypothesis.
  RooRealVar* mu = wspace->var("masterSignal"); 
  RooArgSet* poi = new RooArgSet(*mu); 
  RooArgSet* nullParams = new RooArgSet("nullParams");
  nullParams->addClone(*mu);
  // here we explicitly set the value of the parameters for the null
  nullParams->setRealValue("masterSignal",0); 

  // Step 5, Create a calculator for doing the hypothesis test.
  // because this is a 
  ProfileLikelihoodCalculator plc( *wspace->data("ObservedNumberCountingData"),
				   *wspace->pdf("TopLevelPdf"), *poi, 0.05, nullParams);
	
  wspace->var("tau_0")->Print();
  wspace->var("tau_1")->Print();

  // Step 7, Use the Calculator to get a HypoTestResult
  HypoTestResult* htr = plc.GetHypoTest();
  cout << "-------------------------------------------------" << endl;
  cout << "The p-value for the null is " << htr->NullPValue() << endl;
  cout << "Corresponding to a signifcance of " << htr->Significance() << endl;
  cout << "-------------------------------------------------\n\n" << endl;

  /* observed case should return:
     -------------------------------------------------
     The p-value for the null is 0.0351669
     Corresponding to a signifcance of 1.80975
     -------------------------------------------------
  */


  //////////////////////////////////////////
  // Confidence Interval Stage

  // Step 8, Here we re-use the ProfileLikelihoodCalculator to return a confidence interval.
  // We need to specify what are our parameters of interest
  RooArgSet* paramsOfInterest = nullParams; // they are the same as before in this case
  plc.SetParameters(*paramsOfInterest);
  LikelihoodInterval* lrint = (LikelihoodInterval*) plc.GetInterval();  // that was easy.
  lrint->SetConfidenceLevel(0.95);

  // Step 9c. Get upper and lower limits
  cout << "lower limit on master signal = " <<   lrint->LowerLimit(*mu ) << endl;
  cout << "upper limit on master signal = " <<   lrint->UpperLimit(*mu ) << endl;

  delete lrint;
  delete htr;
  delete wspace;
  delete nullParams;
  delete poi; 

  
}


void rs_numberCountingCombination_observedWithTau()
{

  /////////////////////////////////////////
  // The same example with observed data in a main
  // measurement and an background-only auxiliary 
  // measurement with a factor tau more background
  // than in the main measurement.

  /////////////////////////////////////////
  // The Model building stage
  /////////////////////////////////////////

  // Step 1, define arrays with signal & bkg expectations and background uncertainties
  // We still need the expectation to relate signal in different channels with the master signal
  Double_t s[2] = {20.,10.};           // expected signal
  
  // Step 2, use a RooStats factory to build a PDF for a 
  // number counting combination and add it to the workspace.
  // We need to give the signal expectation to relate the masterSignal
  // to the signal contribution in the individual channels.
  // The model neglects correlations in background uncertainty, 
  // but they could be added without much change to the example.
  NumberCountingPdfFactory f;
  RooWorkspace* wspace = new RooWorkspace();
  f.AddModel(s,2,wspace,"TopLevelPdf", "masterSignal"); 

  // Step 3, use a RooStats factory to add datasets to the workspace.
  // Add the observed data to the workspace in the on-off problem.
  Double_t mainMeas[2] = {123.,117.};      // observed main measurement
  Double_t sideband[2] = {11123.,9876.};    // observed sideband
  Double_t tau[2] = {100.,100.}; // ratio of bkg in sideband to bkg in main measurement, from experimental design.
  f.AddDataWithSideband(mainMeas, sideband, tau, 2, wspace,"ObservedNumberCountingDataWithSideband");

  // see below for a printout of the workspace
  //  wspace->Print();  //uncomment to see structure of workspace

  /////////////////////////////////////////
  // The Hypothesis testing stage:
  /////////////////////////////////////////
  // Step 4, Define the null hypothesis for the calculator
  // Here you need to know the name of the variables corresponding to hypothesis.
  RooRealVar* mu = wspace->var("masterSignal"); 
  RooArgSet* poi = new RooArgSet(*mu); 
  RooArgSet* nullParams = new RooArgSet("nullParams");
  nullParams->addClone(*mu);
  // here we explicitly set the value of the parameters for the null
  nullParams->setRealValue("masterSignal",0); 

  // Step 5, Create a calculator for doing the hypothesis test.
  // because this is a 
  ProfileLikelihoodCalculator plc( *wspace->data("ObservedNumberCountingDataWithSideband"),
				   *wspace->pdf("TopLevelPdf"), *poi, 0.05, nullParams);
				  


  // Step 7, Use the Calculator to get a HypoTestResult
  HypoTestResult* htr = plc.GetHypoTest();
  cout << "-------------------------------------------------" << endl;
  cout << "The p-value for the null is " << htr->NullPValue() << endl;
  cout << "Corresponding to a signifcance of " << htr->Significance() << endl;
  cout << "-------------------------------------------------\n\n" << endl;

  /* observed case should return:
     -------------------------------------------------
     The p-value for the null is 0.0352035
     Corresponding to a signifcance of 1.80928
     -------------------------------------------------
  */


  //////////////////////////////////////////
  // Confidence Interval Stage

  // Step 8, Here we re-use the ProfileLikelihoodCalculator to return a confidence interval.
  // We need to specify what are our parameters of interest
  RooArgSet* paramsOfInterest = nullParams; // they are the same as before in this case
  plc.SetParameters(*paramsOfInterest);
  LikelihoodInterval* lrint = (LikelihoodInterval*) plc.GetInterval();  // that was easy.
  lrint->SetConfidenceLevel(0.95);

  

  // Step 9c. Get upper and lower limits
  cout << "lower limit on master signal = " <<   lrint->LowerLimit(*mu ) << endl;
  cout << "upper limit on master signal = " <<   lrint->UpperLimit(*mu ) << endl;

  delete lrint;
  delete htr;
  delete wspace;
  delete nullParams;
  delete poi; 

  
}
 rs_numberCountingCombination.C:1
 rs_numberCountingCombination.C:2
 rs_numberCountingCombination.C:3
 rs_numberCountingCombination.C:4
 rs_numberCountingCombination.C:5
 rs_numberCountingCombination.C:6
 rs_numberCountingCombination.C:7
 rs_numberCountingCombination.C:8
 rs_numberCountingCombination.C:9
 rs_numberCountingCombination.C:10
 rs_numberCountingCombination.C:11
 rs_numberCountingCombination.C:12
 rs_numberCountingCombination.C:13
 rs_numberCountingCombination.C:14
 rs_numberCountingCombination.C:15
 rs_numberCountingCombination.C:16
 rs_numberCountingCombination.C:17
 rs_numberCountingCombination.C:18
 rs_numberCountingCombination.C:19
 rs_numberCountingCombination.C:20
 rs_numberCountingCombination.C:21
 rs_numberCountingCombination.C:22
 rs_numberCountingCombination.C:23
 rs_numberCountingCombination.C:24
 rs_numberCountingCombination.C:25
 rs_numberCountingCombination.C:26
 rs_numberCountingCombination.C:27
 rs_numberCountingCombination.C:28
 rs_numberCountingCombination.C:29
 rs_numberCountingCombination.C:30
 rs_numberCountingCombination.C:31
 rs_numberCountingCombination.C:32
 rs_numberCountingCombination.C:33
 rs_numberCountingCombination.C:34
 rs_numberCountingCombination.C:35
 rs_numberCountingCombination.C:36
 rs_numberCountingCombination.C:37
 rs_numberCountingCombination.C:38
 rs_numberCountingCombination.C:39
 rs_numberCountingCombination.C:40
 rs_numberCountingCombination.C:41
 rs_numberCountingCombination.C:42
 rs_numberCountingCombination.C:43
 rs_numberCountingCombination.C:44
 rs_numberCountingCombination.C:45
 rs_numberCountingCombination.C:46
 rs_numberCountingCombination.C:47
 rs_numberCountingCombination.C:48
 rs_numberCountingCombination.C:49
 rs_numberCountingCombination.C:50
 rs_numberCountingCombination.C:51
 rs_numberCountingCombination.C:52
 rs_numberCountingCombination.C:53
 rs_numberCountingCombination.C:54
 rs_numberCountingCombination.C:55
 rs_numberCountingCombination.C:56
 rs_numberCountingCombination.C:57
 rs_numberCountingCombination.C:58
 rs_numberCountingCombination.C:59
 rs_numberCountingCombination.C:60
 rs_numberCountingCombination.C:61
 rs_numberCountingCombination.C:62
 rs_numberCountingCombination.C:63
 rs_numberCountingCombination.C:64
 rs_numberCountingCombination.C:65
 rs_numberCountingCombination.C:66
 rs_numberCountingCombination.C:67
 rs_numberCountingCombination.C:68
 rs_numberCountingCombination.C:69
 rs_numberCountingCombination.C:70
 rs_numberCountingCombination.C:71
 rs_numberCountingCombination.C:72
 rs_numberCountingCombination.C:73
 rs_numberCountingCombination.C:74
 rs_numberCountingCombination.C:75
 rs_numberCountingCombination.C:76
 rs_numberCountingCombination.C:77
 rs_numberCountingCombination.C:78
 rs_numberCountingCombination.C:79
 rs_numberCountingCombination.C:80
 rs_numberCountingCombination.C:81
 rs_numberCountingCombination.C:82
 rs_numberCountingCombination.C:83
 rs_numberCountingCombination.C:84
 rs_numberCountingCombination.C:85
 rs_numberCountingCombination.C:86
 rs_numberCountingCombination.C:87
 rs_numberCountingCombination.C:88
 rs_numberCountingCombination.C:89
 rs_numberCountingCombination.C:90
 rs_numberCountingCombination.C:91
 rs_numberCountingCombination.C:92
 rs_numberCountingCombination.C:93
 rs_numberCountingCombination.C:94
 rs_numberCountingCombination.C:95
 rs_numberCountingCombination.C:96
 rs_numberCountingCombination.C:97
 rs_numberCountingCombination.C:98
 rs_numberCountingCombination.C:99
 rs_numberCountingCombination.C:100
 rs_numberCountingCombination.C:101
 rs_numberCountingCombination.C:102
 rs_numberCountingCombination.C:103
 rs_numberCountingCombination.C:104
 rs_numberCountingCombination.C:105
 rs_numberCountingCombination.C:106
 rs_numberCountingCombination.C:107
 rs_numberCountingCombination.C:108
 rs_numberCountingCombination.C:109
 rs_numberCountingCombination.C:110
 rs_numberCountingCombination.C:111
 rs_numberCountingCombination.C:112
 rs_numberCountingCombination.C:113
 rs_numberCountingCombination.C:114
 rs_numberCountingCombination.C:115
 rs_numberCountingCombination.C:116
 rs_numberCountingCombination.C:117
 rs_numberCountingCombination.C:118
 rs_numberCountingCombination.C:119
 rs_numberCountingCombination.C:120
 rs_numberCountingCombination.C:121
 rs_numberCountingCombination.C:122
 rs_numberCountingCombination.C:123
 rs_numberCountingCombination.C:124
 rs_numberCountingCombination.C:125
 rs_numberCountingCombination.C:126
 rs_numberCountingCombination.C:127
 rs_numberCountingCombination.C:128
 rs_numberCountingCombination.C:129
 rs_numberCountingCombination.C:130
 rs_numberCountingCombination.C:131
 rs_numberCountingCombination.C:132
 rs_numberCountingCombination.C:133
 rs_numberCountingCombination.C:134
 rs_numberCountingCombination.C:135
 rs_numberCountingCombination.C:136
 rs_numberCountingCombination.C:137
 rs_numberCountingCombination.C:138
 rs_numberCountingCombination.C:139
 rs_numberCountingCombination.C:140
 rs_numberCountingCombination.C:141
 rs_numberCountingCombination.C:142
 rs_numberCountingCombination.C:143
 rs_numberCountingCombination.C:144
 rs_numberCountingCombination.C:145
 rs_numberCountingCombination.C:146
 rs_numberCountingCombination.C:147
 rs_numberCountingCombination.C:148
 rs_numberCountingCombination.C:149
 rs_numberCountingCombination.C:150
 rs_numberCountingCombination.C:151
 rs_numberCountingCombination.C:152
 rs_numberCountingCombination.C:153
 rs_numberCountingCombination.C:154
 rs_numberCountingCombination.C:155
 rs_numberCountingCombination.C:156
 rs_numberCountingCombination.C:157
 rs_numberCountingCombination.C:158
 rs_numberCountingCombination.C:159
 rs_numberCountingCombination.C:160
 rs_numberCountingCombination.C:161
 rs_numberCountingCombination.C:162
 rs_numberCountingCombination.C:163
 rs_numberCountingCombination.C:164
 rs_numberCountingCombination.C:165
 rs_numberCountingCombination.C:166
 rs_numberCountingCombination.C:167
 rs_numberCountingCombination.C:168
 rs_numberCountingCombination.C:169
 rs_numberCountingCombination.C:170
 rs_numberCountingCombination.C:171
 rs_numberCountingCombination.C:172
 rs_numberCountingCombination.C:173
 rs_numberCountingCombination.C:174
 rs_numberCountingCombination.C:175
 rs_numberCountingCombination.C:176
 rs_numberCountingCombination.C:177
 rs_numberCountingCombination.C:178
 rs_numberCountingCombination.C:179
 rs_numberCountingCombination.C:180
 rs_numberCountingCombination.C:181
 rs_numberCountingCombination.C:182
 rs_numberCountingCombination.C:183
 rs_numberCountingCombination.C:184
 rs_numberCountingCombination.C:185
 rs_numberCountingCombination.C:186
 rs_numberCountingCombination.C:187
 rs_numberCountingCombination.C:188
 rs_numberCountingCombination.C:189
 rs_numberCountingCombination.C:190
 rs_numberCountingCombination.C:191
 rs_numberCountingCombination.C:192
 rs_numberCountingCombination.C:193
 rs_numberCountingCombination.C:194
 rs_numberCountingCombination.C:195
 rs_numberCountingCombination.C:196
 rs_numberCountingCombination.C:197
 rs_numberCountingCombination.C:198
 rs_numberCountingCombination.C:199
 rs_numberCountingCombination.C:200
 rs_numberCountingCombination.C:201
 rs_numberCountingCombination.C:202
 rs_numberCountingCombination.C:203
 rs_numberCountingCombination.C:204
 rs_numberCountingCombination.C:205
 rs_numberCountingCombination.C:206
 rs_numberCountingCombination.C:207
 rs_numberCountingCombination.C:208
 rs_numberCountingCombination.C:209
 rs_numberCountingCombination.C:210
 rs_numberCountingCombination.C:211
 rs_numberCountingCombination.C:212
 rs_numberCountingCombination.C:213
 rs_numberCountingCombination.C:214
 rs_numberCountingCombination.C:215
 rs_numberCountingCombination.C:216
 rs_numberCountingCombination.C:217
 rs_numberCountingCombination.C:218
 rs_numberCountingCombination.C:219
 rs_numberCountingCombination.C:220
 rs_numberCountingCombination.C:221
 rs_numberCountingCombination.C:222
 rs_numberCountingCombination.C:223
 rs_numberCountingCombination.C:224
 rs_numberCountingCombination.C:225
 rs_numberCountingCombination.C:226
 rs_numberCountingCombination.C:227
 rs_numberCountingCombination.C:228
 rs_numberCountingCombination.C:229
 rs_numberCountingCombination.C:230
 rs_numberCountingCombination.C:231
 rs_numberCountingCombination.C:232
 rs_numberCountingCombination.C:233
 rs_numberCountingCombination.C:234
 rs_numberCountingCombination.C:235
 rs_numberCountingCombination.C:236
 rs_numberCountingCombination.C:237
 rs_numberCountingCombination.C:238
 rs_numberCountingCombination.C:239
 rs_numberCountingCombination.C:240
 rs_numberCountingCombination.C:241
 rs_numberCountingCombination.C:242
 rs_numberCountingCombination.C:243
 rs_numberCountingCombination.C:244
 rs_numberCountingCombination.C:245
 rs_numberCountingCombination.C:246
 rs_numberCountingCombination.C:247
 rs_numberCountingCombination.C:248
 rs_numberCountingCombination.C:249
 rs_numberCountingCombination.C:250
 rs_numberCountingCombination.C:251
 rs_numberCountingCombination.C:252
 rs_numberCountingCombination.C:253
 rs_numberCountingCombination.C:254
 rs_numberCountingCombination.C:255
 rs_numberCountingCombination.C:256
 rs_numberCountingCombination.C:257
 rs_numberCountingCombination.C:258
 rs_numberCountingCombination.C:259
 rs_numberCountingCombination.C:260
 rs_numberCountingCombination.C:261
 rs_numberCountingCombination.C:262
 rs_numberCountingCombination.C:263
 rs_numberCountingCombination.C:264
 rs_numberCountingCombination.C:265
 rs_numberCountingCombination.C:266
 rs_numberCountingCombination.C:267
 rs_numberCountingCombination.C:268
 rs_numberCountingCombination.C:269
 rs_numberCountingCombination.C:270
 rs_numberCountingCombination.C:271
 rs_numberCountingCombination.C:272
 rs_numberCountingCombination.C:273
 rs_numberCountingCombination.C:274
 rs_numberCountingCombination.C:275
 rs_numberCountingCombination.C:276
 rs_numberCountingCombination.C:277
 rs_numberCountingCombination.C:278
 rs_numberCountingCombination.C:279
 rs_numberCountingCombination.C:280
 rs_numberCountingCombination.C:281
 rs_numberCountingCombination.C:282
 rs_numberCountingCombination.C:283
 rs_numberCountingCombination.C:284
 rs_numberCountingCombination.C:285
 rs_numberCountingCombination.C:286
 rs_numberCountingCombination.C:287
 rs_numberCountingCombination.C:288
 rs_numberCountingCombination.C:289
 rs_numberCountingCombination.C:290
 rs_numberCountingCombination.C:291
 rs_numberCountingCombination.C:292
 rs_numberCountingCombination.C:293
 rs_numberCountingCombination.C:294
 rs_numberCountingCombination.C:295
 rs_numberCountingCombination.C:296
 rs_numberCountingCombination.C:297
 rs_numberCountingCombination.C:298
 rs_numberCountingCombination.C:299
 rs_numberCountingCombination.C:300
 rs_numberCountingCombination.C:301
 rs_numberCountingCombination.C:302
 rs_numberCountingCombination.C:303
 rs_numberCountingCombination.C:304
 rs_numberCountingCombination.C:305
 rs_numberCountingCombination.C:306
 rs_numberCountingCombination.C:307
 rs_numberCountingCombination.C:308
 rs_numberCountingCombination.C:309
 rs_numberCountingCombination.C:310
 rs_numberCountingCombination.C:311
 rs_numberCountingCombination.C:312
 rs_numberCountingCombination.C:313
 rs_numberCountingCombination.C:314
 rs_numberCountingCombination.C:315
 rs_numberCountingCombination.C:316
 rs_numberCountingCombination.C:317
 rs_numberCountingCombination.C:318
 rs_numberCountingCombination.C:319
 rs_numberCountingCombination.C:320
 rs_numberCountingCombination.C:321
 rs_numberCountingCombination.C:322
 rs_numberCountingCombination.C:323
 rs_numberCountingCombination.C:324
 rs_numberCountingCombination.C:325
 rs_numberCountingCombination.C:326
 rs_numberCountingCombination.C:327
 rs_numberCountingCombination.C:328
 rs_numberCountingCombination.C:329
 rs_numberCountingCombination.C:330
 rs_numberCountingCombination.C:331
 rs_numberCountingCombination.C:332
 rs_numberCountingCombination.C:333
 rs_numberCountingCombination.C:334
 rs_numberCountingCombination.C:335
 rs_numberCountingCombination.C:336
 rs_numberCountingCombination.C:337
 rs_numberCountingCombination.C:338
 rs_numberCountingCombination.C:339
 rs_numberCountingCombination.C:340
 rs_numberCountingCombination.C:341
 rs_numberCountingCombination.C:342
 rs_numberCountingCombination.C:343
 rs_numberCountingCombination.C:344
 rs_numberCountingCombination.C:345
 rs_numberCountingCombination.C:346
 rs_numberCountingCombination.C:347
 rs_numberCountingCombination.C:348
 rs_numberCountingCombination.C:349
 rs_numberCountingCombination.C:350
 rs_numberCountingCombination.C:351
 rs_numberCountingCombination.C:352
 rs_numberCountingCombination.C:353
 rs_numberCountingCombination.C:354
 rs_numberCountingCombination.C:355
 rs_numberCountingCombination.C:356
 rs_numberCountingCombination.C:357
 rs_numberCountingCombination.C:358
 rs_numberCountingCombination.C:359
 rs_numberCountingCombination.C:360
 rs_numberCountingCombination.C:361
 rs_numberCountingCombination.C:362
 rs_numberCountingCombination.C:363
 rs_numberCountingCombination.C:364
 rs_numberCountingCombination.C:365
 rs_numberCountingCombination.C:366
 rs_numberCountingCombination.C:367
 rs_numberCountingCombination.C:368
 rs_numberCountingCombination.C:369
 rs_numberCountingCombination.C:370
 rs_numberCountingCombination.C:371
 rs_numberCountingCombination.C:372
 rs_numberCountingCombination.C:373
 rs_numberCountingCombination.C:374
 rs_numberCountingCombination.C:375
 rs_numberCountingCombination.C:376
 rs_numberCountingCombination.C:377
 rs_numberCountingCombination.C:378
 rs_numberCountingCombination.C:379
 rs_numberCountingCombination.C:380
 rs_numberCountingCombination.C:381
 rs_numberCountingCombination.C:382
 rs_numberCountingCombination.C:383
 rs_numberCountingCombination.C:384
 rs_numberCountingCombination.C:385
 rs_numberCountingCombination.C:386
 rs_numberCountingCombination.C:387
 rs_numberCountingCombination.C:388
 rs_numberCountingCombination.C:389
 rs_numberCountingCombination.C:390
 rs_numberCountingCombination.C:391
 rs_numberCountingCombination.C:392
 rs_numberCountingCombination.C:393
 rs_numberCountingCombination.C:394
 rs_numberCountingCombination.C:395
 rs_numberCountingCombination.C:396
 rs_numberCountingCombination.C:397
 rs_numberCountingCombination.C:398
 rs_numberCountingCombination.C:399
 rs_numberCountingCombination.C:400
 rs_numberCountingCombination.C:401
 rs_numberCountingCombination.C:402
 rs_numberCountingCombination.C:403
 rs_numberCountingCombination.C:404
 rs_numberCountingCombination.C:405
 rs_numberCountingCombination.C:406
 rs_numberCountingCombination.C:407
 rs_numberCountingCombination.C:408
 rs_numberCountingCombination.C:409
 rs_numberCountingCombination.C:410
 rs_numberCountingCombination.C:411
 rs_numberCountingCombination.C:412
 rs_numberCountingCombination.C:413
 rs_numberCountingCombination.C:414
 rs_numberCountingCombination.C:415
 rs_numberCountingCombination.C:416
 rs_numberCountingCombination.C:417
 rs_numberCountingCombination.C:418
 rs_numberCountingCombination.C:419