Merge pull request #3 from jingyuluo/master

Add stability scripts

Studies/Draw_PCCratio.C

@@ -0,0 +1,107 @@
+void Draw_PCCratio(){
+
+    gStyle->SetOptStat(0);
+
+    TCanvas *c = new TCanvas("c", "c", 800, 500);
+    c->cd();
+    c->SetTickx();
+    c->SetTicky();
+
+
+    TFile *f = new TFile("Overall.root");
+
+    TH1F *h_PCC_layer_0 = (TH1F*)f->Get("h_PCC_Layer_1");
+    TH1F *h_PCC_layer_1 = (TH1F*)f->Get("h_PCC_Layer_2");
+    TH1F *h_PCC_layer_2 = (TH1F*)f->Get("h_PCC_Layer_3");
+    TH1F *h_PCC_layer_3 = (TH1F*)f->Get("h_PCC_Layer_4");
+    TH1F *h_PCC_total = (TH1F*)f->Get("h_PCC_total");
+
+    TH1F *ratio_0 = (TH1F*)h_PCC_layer_0->Clone();
+    TH1F *ratio_1 = (TH1F*)h_PCC_layer_1->Clone();
+    TH1F *ratio_2 = (TH1F*)h_PCC_layer_2->Clone();
+    TH1F *ratio_3 = (TH1F*)h_PCC_layer_3->Clone();
+
+    ratio_0->Divide(h_PCC_total);
+    ratio_1->Divide(h_PCC_total);
+    ratio_2->Divide(h_PCC_total);
+    ratio_3->Divide(h_PCC_total);
+
+    for(int i=0; i<ratio_0->GetNbinsX(); i++){
+        ratio_0->SetBinError(i, 0.00001);
+
+    }
+
+    for(int j=0; j<ratio_1->GetNbinsX(); j++){
+        ratio_1->SetBinError(j, 0.00001);
+
+    }
+
+    for(int k=0; k<ratio_2->GetNbinsX(); k++){
+        ratio_2->SetBinError(k, 0.00001);
+
+    } 
+
+    for(int l=0; l<ratio_3->GetNbinsX(); l++){
+        ratio_3->SetBinError(l, 0.00001);
+
+    }
+
+    ratio_0->SetTitle("");
+
+    ratio_0->SetMarkerStyle(23);
+    ratio_1->SetMarkerStyle(23);
+    ratio_2->SetMarkerStyle(23);
+    ratio_3->SetMarkerStyle(23);
+
+    ratio_0->SetMarkerColor(kBlue);
+    ratio_1->SetMarkerColor(kRed);
+    ratio_2->SetMarkerColor(kGreen);
+    ratio_3->SetMarkerColor(6);
+
+    ratio_0->SetLineColor(kBlue);
+    ratio_1->SetLineColor(kRed);
+    ratio_2->SetLineColor(kGreen);
+    ratio_3->SetLineColor(6);
+
+  //  ratio_0->SetLineWidth(4);
+  //  ratio_1->SetLineWidth(4);
+  //  ratio_2->SetLineWidth(4);
+  //  ratio_3->SetLineWidth(4);
+
+    ratio_0->GetXaxis()->SetTitle("Date (Day/Month)");
+    ratio_0->GetXaxis()->SetTimeDisplay(1);
+    //ratio_0->GetXaxis()->SetTimeOffset("GMT");
+    ratio_0->GetXaxis()->SetNdivisions(-503);
+    ratio_0->GetXaxis()->SetTimeFormat("  %d/%m %F1970-01-01 00:00:00");
+    ratio_0->GetYaxis()->SetTitle("Relative Contribution");
+    //ratio_0->GetXaxis()->SetRangeUser(251000, 261000);
+    ratio_0->GetYaxis()->SetRangeUser(0.17,0.35);
+    ratio_0->Draw("PE");
+    ratio_1->Draw("PE SAME");
+    ratio_2->Draw("PE SAME");
+    ratio_3->Draw("PE SAME");
+
+    TLegend *len = new TLegend(0.5, 0.4, 0.85, 0.55);
+    len->SetFillColor(0);
+    len->SetLineColor(0);
+    len->AddEntry(ratio_0,"barrel layer 2", "P");
+    len->AddEntry(ratio_1,"barrel layer 3", "P");
+    len->AddEntry(ratio_2,"forward disk 1", "P");
+    len->AddEntry(ratio_3,"forward disk 2", "P");
+
+
+    len->Draw("SAME");
+
+    TLatex *text=new TLatex(0.72,0.85,"2015  (13TeV)");
+    text->SetNDC();
+    text->SetTextFont(62);
+    text->SetTextSize(0.05);
+    TLatex *text2=new TLatex(0.15,0.85,"CMS #bf{#scale[0.75]{#it{Preliminary}}}");
+    text2->SetNDC();
+    text2->SetTextSize(0.05);
+    text2->SetTextFont(62);
+    text->Draw("same");
+    text2->Draw("same");
+    c->SaveAs("new_PCC_stability.eps");
+    c->SaveAs("new_PCC_stability.png");
+}

Studies/PCCratioJobs.py

@@ -0,0 +1,75 @@
+import ROOT
+import sys
+import argparse
+import os
+import subprocess
+
+parser = argparse.ArgumentParser()
+parser.add_argument("-p", "--path", help="EOS path to certtree files")
+parser.add_argument("-d", "--dir", default="JobDir", help="Output Directory")
+parser.add_argument("-o", "--outPath", default="", help="Specify the path of output files")
+parser.add_argument("-j", "--json", default="", help="JSON file fo run/LSs to filter")
+parser.add_argument("-s", "--sub", default=False, action="store_true", help="bsub created jobs")
+
+args=parser.parse_args()
+
+def MakeJob(outputdir, jobid, filename):
+
+    joblines=[]
+    joblines.append("source /cvmfs/cms.cern.ch/cmsset_default.sh")
+    joblines.append("cd "+outputdir)
+    joblines.append("cmsenv")
+    makeDataCMD ="python ../PCCstability.py --certfile="+args.path+"/"+filename
+    if args.json!="":
+        makeDataCMD=makeDataCMD+" --json="+args.json
+
+    if args.outPath!="":
+        makeDataCMD=makeDataCMD+" --output="+args.outPath
+
+    makeDataCMD=makeDataCMD+" --label="+str(jobid)
+
+    joblines.append(makeDataCMD)
+
+    scriptFile=open(outputdir+"/job_"+str(jobid)+".sh", "w+")
+    for line in joblines:
+        scriptFile.write(line+"\n")
+
+    scriptFile.close()
+
+def SubmitJob(job, queue="1nh"):
+    baseName=str(job.split(".")[0])
+    cmd="bsub -q "+queue+" -J "+baseName+" -o "+baseName+".log < "+str(job)
+    output=os.system(cmd)
+    if output!=0:
+        print job, "did not submit properly"
+        print cmd
+
+
+fileinfos=subprocess.check_output(["/afs/cern.ch/project/eos/installation/0.3.4/bin/eos.select","ls", args.path])
+fileinfos=fileinfos.split("\n")
+
+filenames={}
+for fileinfo in fileinfos:
+
+    filename=fileinfo
+    if filename.find(".root") == -1:
+        continue
+    jobid=filename.split("/")[-1].split(".")[0].split("_")[-1]
+
+    filenames[int(jobid)]=filename
+
+fullOutPath = os.getcwd()
+if not os.path.exists(args.dir):
+
+    os.makedirs(args.dir)
+fullOutPath=fullOutPath+"/"+args.dir
+
+for job in filenames:
+    MakeJob(fullOutPath, job, filenames[job])
+
+if args.sub: 
+    print "Submitting", len(filename), "jobs"
+
+    for job in filenames:
+
+        SubmitJob(args.dir+"/job_"+str(job)+".sh", "1nh")

Studies/PCCstability.py

@@ -0,0 +1,93 @@
+import ROOT
+from ROOT import TH1F
+import sys,os
+import numpy
+import math
+import argparse
+import pickle
+import time
+import json
+
+parser = argparse.ArgumentParser(description='Process some integers')
+parser.add_argument('-c','--certfile', type=str, default="", help='The certfile to input')
+parser.add_argument('-o','--output',type=str, default="", help='The path for the output file')
+parser.add_argument('-j','--json', type=str, default="", help='JSON file of run/LSs to filter')
+parser.add_argument('-l','--label', type=str, default="", help='The label of output file')
+
+args=parser.parse_args()
+
+if args.output!="":
+    outpath=args.output
+
+    if outpath.find("/store")==0:
+        outpath="root://eoscms//eos/cms"+outpath
+    newfilename="PCC_ratio_"+args.label+".root"
+    newfile=ROOT.TFile.Open(newfilename, "recreate")
+
+
+#hist_PCCratio_Layer = {}
+hist_PCC_layer = {}
+
+for iLayer in range(5):
+    #hist_PCCratio_Layer[iLayer] = TH1F("h_PCCratio_Layer_"+str(iLayer), "h_PCCratio_Layer_"+str(iLayer), 10000, 0, 10000)
+    #hist_PCCratio_Layer[iLayer].Sumw2()
+
+    hist_PCC_layer[iLayer] = TH1F("h_PCC_Layer_"+str(iLayer), "h_PCC_Layer_"+str(iLayer), 800, 1436000000, 1437000000)
+    #hist_PCC_layer[iLayer].Sumw2()
+hist_PCC_total = TH1F("h_PCC_total", "h_PCC_total", 800, 1436000000, 1437000000)
+hist_PCC_nLS = TH1F("h_PCC_nLS", "h_PCC_nLS", 800, 1436000000, 1437000000)
+
+runLSData = json.load(open(args.json))
+
+
+def ISRunLSInList(run, LS):
+    if not runLSData.has_key(str(run)):
+        return False
+
+    else:
+        for LSRange in runLSData[str(run)]:
+            if int(LSRange[0])<=LS and int(LSRange[1])>=LS:
+                return True 
+
+if args.certfile!="":
+    filename=args.certfile
+
+    if filename.find("/store")==0:
+        filename="root://eoscms//eos/cms"+filename
+
+    tfile=ROOT.TFile.Open(filename)
+
+    tree=tfile.Get("certtree")
+    tree.SetBranchStatus("*", 0)
+    tree.SetBranchStatus("run*", 1)
+    tree.SetBranchStatus("LS*", 1)
+    tree.SetBranchStatus("timeStamp*", 1)
+    tree.SetBranchStatus("nCluster*", 1)
+    tree.SetBranchStatus("nPCPerLayer*", 1)
+
+    nentries = tree.GetEntries()
+
+    for iev in range(nentries):
+
+        tree.GetEntry(iev)
+        if iev%1==0:
+            print "iev,", iev
+
+        if ISRunLSInList(tree.run, tree.LS):
+            print "yes!"
+            print tree.timeStamp
+           # if hist_PCC_nLS.GetBinContent(hist_PCC_nLS.FindBin(tree.run))<50:
+            hist_PCC_total.Fill(tree.timeStamp, tree.nCluster)
+            hist_PCC_nLS.Fill(tree.timeStamp, 1)
+
+            for iLayer in range(5):
+                hist_PCC_layer[iLayer].Fill(tree.timeStamp, tree.nPCPerLayer[iLayer])    
+
+if args.output!="":
+    #newfile.WriteTObject(
+    newfile.Write()
+    newfile.Close()
+
+
+
+