Muon 2nd control region

config/cross_section_config.py

@@ -356,7 +356,7 @@ def __fill_defaults__( self ):
 
         self.tree_path_control_templates = {
                                     'electron' : 'TTbar_plus_X_analysis/EPlusJets/QCD non iso e+jets/FitVariables',
-                                    'muon' : 'TTbar_plus_X_analysis/MuPlusJets/QCD non iso mu+jets/FitVariables'
+                                    'muon' : 'TTbar_plus_X_analysis/MuPlusJets/QCD iso > 0.3/FitVariables'
                                     }
         self.variable_path_templates = {
                             'MET' : 'TTbar_plus_X_analysis/{channel}/{selection}/FitVariables/MET',
@@ -379,8 +379,8 @@ def __fill_defaults__( self ):
         self.muon_control_region = 'QCD non iso mu+jets ge3j'
         self.muon_control_region_systematic = 'QCD non iso mu+jets ge3j'  # no systematic yet
         if self.centre_of_mass_energy == 13:
-            self.muon_control_region = 'QCD non iso mu+jets'
-            self.muon_control_region_systematic = 'QCD non iso mu+jets'
+            self.muon_control_region = 'QCD iso > 0.3'
+            self.muon_control_region_systematic = 'QCD 0.12 < iso <= 0.3'
 
         self.include_higgs = False
 

experimental/add_control_region.py

@@ -0,0 +1,69 @@
+'''
+Created on 26 Aug 2015
+
+@author: kreczko
+'''
+from rootpy.io.file import root_open, File
+from tools.ROOT_utils import root_mkdir
+import shutil
+import subprocess
+
+input_folder = '/hdfs/TopQuarkGroup/run2/atOutput/13TeV/50ns/'
+output_folder = input_folder.replace('/hdfs', '')
+input_files = ['data_electron_tree.root', 'data_muon_tree.root', 
+               'SingleTop_tree.root', 'TTJets_PowhegPythia8_tree.root', 'VJets_tree.root', 'QCD_Muon_tree.root']
+
+
+def create_new_trees(input_file, suffix = ''):
+    tree1_name = 'TTbar_plus_X_analysis/MuPlusJets/QCD non iso mu+jets/FitVariables' + suffix
+    tree2_name = 'TTbar_plus_X_analysis/MuPlusJets/QCD non iso mu+jets/Muon/Muons' + suffix
+    s_cr1 = 'relIso_04_deltaBeta <= 0.3 && relIso_04_deltaBeta > 0.12'
+    s_cr2 = 'relIso_04_deltaBeta > 0.3'
+    cr1 = 'TTbar_plus_X_analysis/MuPlusJets/QCD 0.12 < iso <= 0.3'
+    cr2 = 'TTbar_plus_X_analysis/MuPlusJets/QCD iso > 0.3'
+
+    with root_open(input_file) as f:
+        t1 = f.Get(tree1_name)
+        t2 = f.Get(tree2_name)
+        t1.AddFriend(t2)
+
+    #     h1 = t1.Draw('MET', 'relIso_04_deltaBeta > 0.3')
+    #     h2 = t1.Draw(
+    #         'MET', 'relIso_04_deltaBeta <= 0.3 && relIso_04_deltaBeta > 0.12')
+    #     h3 = t1.Draw('MET', 'relIso_04_deltaBeta > 0.12')
+    #     h4 = t2.Draw('relIso_04_deltaBeta', 'relIso_04_deltaBeta > 0.12')
+    #     print h1.integral()
+    #     print h2.integral()
+    #     print h3.integral(), h1.integral() + h2.integral()
+
+        output = File('test.root', 'recreate')
+        output.mkdir(cr1, recurse=True)
+        output.mkdir(cr2, recurse=True)
+        output.cd(cr2)
+        new_tree1 = t1.CopyTree(s_cr2)
+        new_tree1.Write()
+        output.cd(cr1)
+        new_tree2 = t1.CopyTree(s_cr1)
+        new_tree2.Write()
+        output.close()
+
+    new_tree1 = None
+    new_tree2 = None
+
+    f_out = File(input_file, 'update')
+    root_mkdir(f_out, cr1)
+    root_mkdir(f_out, cr2)
+    with root_open('test.root') as f_in:
+        f_out.cd(cr1)
+        new_tree1 = f_in.Get(cr1 + '/FitVariables' + suffix).CloneTree()
+        f_out.cd(cr2)
+        new_tree2 = f_in.Get(cr2 + '/FitVariables' + suffix).CloneTree()
+
+for f in input_files:
+    shutil.copy(input_folder + f, f)
+    for suffix in ['', '_JERDown', '_JERUp', '_JESDown', '_JESUp']:
+        if 'data' in f and not suffix == '':
+            continue
+        create_new_trees(f, suffix=suffix)
+    subprocess.call(['hadoop', 'fs', '-rm', output_folder + f])
+    subprocess.call(['hadoop', 'fs', '-copyFromLocal', f, output_folder + f])

src/cross_section_measurement/create_measurement.py

@@ -221,17 +221,25 @@ def create_measurement(com, category, variable, channel, phase_space, norm_metho
     m.addShapeForSample('QCD', m_qcd, False)
     norm_qcd = deepcopy(m_qcd)
     # we want QCD shape and normalisation to be separate
-    if category == 'QCD_shape' and channel == 'electron':
+    if category == 'QCD_shape':
         for sample in norm_qcd.samples.keys():
             tree = norm_qcd.samples[sample]['input'].tree_name
-            tree = tree.replace(config.electron_control_region_systematic,
-                                config.electron_control_region)
+            if channel == 'electron':
+                tree = tree.replace(config.electron_control_region_systematic,
+                                    config.electron_control_region)
+            else:
+                tree = tree.replace(config.muon_control_region_systematic,
+                                    config.muon_control_region)
             norm_qcd.samples[sample]['input'].tree_name = tree
-    if 'QCD_cross_section' in category and channel == 'electron':
+    if 'QCD_cross_section' in category:
         for sample in norm_qcd.samples.keys():
             tree = norm_qcd.samples[sample]['input'].tree_name
-            tree = tree.replace(config.electron_control_region,
-                                config.electron_control_region_systematic)
+            if channel == 'electron':
+                tree = tree.replace(config.electron_control_region,
+                                    config.electron_control_region_systematic)
+            else:
+                tree = tree.replace(config.muon_control_region,
+                                    config.muon_control_region_systematic)
             norm_qcd.samples[sample]['input'].tree_name = tree
 
     m.addNormForSample('QCD', norm_qcd, False)

src/cross_section_measurement/make_control_plots_fromTrees.py

@@ -83,7 +83,7 @@ def make_plot( channel, x_axis_title, y_axis_title,
         if normalise_to_fit:
             normalisation = normalisations_muon[norm_variable]
         if use_qcd_data_region:
-            qcd_data_region = 'QCD non iso mu+jets'
+            qcd_data_region = 'QCD iso > 0.3'
         if not 'QCD' in channel and not 'NPU' in branchName:
             weightBranchSignalRegion += ' * MuonEfficiencyCorrection'
 
@@ -808,7 +808,8 @@ def make_plot( channel, x_axis_title, y_axis_title,
     for channel, label in {
                             'electronQCDNonIso' : 'EPlusJets/QCD non iso e+jets',
                             'electronQCDConversions' : 'EPlusJets/QCDConversions', 
-                            'muonQCDNonIso' : 'MuPlusJets/QCD non iso mu+jets'
+                            'muonQCDNonIso' : 'MuPlusJets/QCD iso > 0.3',
+                            'muonQCDNonIso2' : 'MuPlusJets/QCD 0.12 < iso <= 0.3',
                             }.iteritems() :
         b_tag_bin = '0btag'
 
@@ -823,7 +824,10 @@ def make_plot( channel, x_axis_title, y_axis_title,
         if channel == 'electronQCDConversions':
             treeName = 'EPlusJets/QCDConversions'
         elif channel == 'muonQCDNonIso':
-            treeName = 'MuPlusJets/QCD non iso mu+jets'
+            treeName = 'MuPlusJets/QCD iso > 0.3'
+            signalTreeName = 'MuPlusJets/Ref selection'
+        elif channel == 'muonQCDNonIso2':
+            treeName = 'MuPlusJets/QCD 0.12 < iso <= 0.3'
             signalTreeName = 'MuPlusJets/Ref selection'
 
         ###################################################

tools/measurement.py

@@ -257,6 +257,7 @@ def read_norm(self, sample):
             self.aux_info_norms[sample]['norm_factor'] = round(ratio, 2)
             self.aux_info_norms[sample]['n_mc_control'] = n_mc_control
             self.aux_info_norms[sample]['n_mc_signal_region'] = n_mc_signal_region
+            self.aux_info_norms[sample]['n_data_control'] = n_data_control
         else:
             meas_log.warning(
                 'No MC entry found for sample "{0}", using control region normalisation'.format(sample))