ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/cbrown/Development/Plotting/Modules/SUSYScan.C

Comparing UserCode/cbrown/Development/Plotting/Modules/SUSYScan.C (file contents):
Revision 1.3 by buchmann, Thu Mar 15 20:34:15 2012 UTC vs.
Revision 1.5 by buchmann, Sun Mar 18 17:53:29 2012 UTC

#	Line 94 | Line 94 | bool initialized_t2=false;
94
95		int srmcpretries=0;
96
97	<	void load_scan_sample(int a, int b, int &scanfileindex,int scantype,bool isretry=false) {
98	<	/*
97	>	void load_scan_sample(int a, int b, int &scanfileindex,int scantype,int mglu, int mlsp, bool isretry=false) {
98	>
99		// There is no need to define your sample here. That is done in Setup.C where you define the loading directory!
100
101	<	dout << "Going to load file with Xzone=" << a << " and  Yzone=" << b << endl;
101	>	dout << "Going to load file with Xzone=" << a << " and  Yzone=" << b << " for scantype (mSugra: " << (bool)(scantype==mSUGRA) << " , SMS: " << (bool)(scantype==SMS) << " , GMSB: " << (bool)(scantype==GMSB) << ")" << endl;
102		stringstream filetoload;
103		char hostname[1023];
104		gethostname(hostname,1023);
105	<
105	>	string samplename;
106		if((Contains(hostname,"t3ui")||Contains(hostname,"t3wn"))) { // local case
107		filetoload << "/shome/buchmann/ntuples/"<<PlottingSetup::ScanSampleDirectory;
108	<	if(scantype==mSUGRA) filetoload << "/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
109	<	if(scantype==SMS) filetoload << "/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
110	<	if(scantype==GMSB) filetoload << "/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
108	>	if(scantype==mSUGRA) {
109	>	//filetoload << "/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
110	>	filetoload.str("");
111	>	filetoload << "/shome/lbaeni/jzb/" << PlottingSetup::ScanSampleDirectory  << "/mSUGRA_M0_" << mglu << "__M12_" << mlsp << ".root";
112	>	samplename="mSUGRA";
113	>	//      filetoload << "/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
114	>	}
115	>	if(scantype==SMS) {
116	>	filetoload << "/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
117	>	samplename="SMS";
118	>	}
119	>	if(scantype==GMSB) {
120	>	filetoload << "/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
121	>	samplename="GMSB";
122	>	}
123		if(scansample.collection.size()<1||!Contains(((scansample.collection)[(scansample.collection).size()-1]).filename,"_"+any2string(a)+"_"+any2string(b))) {
124		dout << "The last sample is NOT the same one as the current one, possibly popping off last one and adding the new one." << endl;
125	<	if((scansample.collection).size()>1) {
125	>	if((scansample.collection).size()>=1) {
126		scansample.RemoveLastSample();
127		}
128		scanfileindex=(scansample.collection).size();
129		//New: Loading file when necessary, not before (avoiding high memory usage and startup times)
130	<	if(scanfileindex!=0) scansample.AddSample(filetoload.str(),"scansample",1,1,false,true,scanfileindex,kRed);
130	>
131	>	scansample.AddSample(filetoload.str(),samplename,1,1,false,true,scanfileindex,kRed);
132	>	dout << "Just added the following file: " << filetoload.str() << endl;
133		} else {
134		dout << "Last sample is the same as the current one. Recycling it." << endl;
135		scanfileindex=(scansample.collection).size()-1;
#	Line 123 | Line 137 | void load_scan_sample(int a, int b, int
137		dout << " Going to use the following file: " << filetoload.str() << endl;
138		} // end of t3 case
139		else {
140	+	//CRAB case
141		write_info(__FUNCTION__,"Hello, CRAB! Might need to load files ...");
142		PlottingSetup::limitpatience=PlottingSetup::limitpatienceCRAB;
143		stringstream copyfile;
144	+	gSystem->Exec("mkdir -p local_storage");
145		copyfile << "lcg-cp -b -T srmv2 srm://storage01.lcg.cscs.ch:8443/srm/managerv2?SFN=/pnfs/lcg.cscs.ch/cms/trivcat/store/user/buchmann/bussola/" << PlottingSetup::ScanSampleDirectory;
146	<	if(scantype==mSUGRA) copyfile << "/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root ";
147	<	if(scantype==SMS) copyfile << "/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root ";
148	<	if(scantype==GMSB) copyfile << "/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root ";
146	>	if(scantype==mSUGRA) {
147	>	copyfile<< "/mSUGRA_M0_" << mglu << "__M12_" << mlsp << ".root local_storage/mSUGRA_M0_" << mglu << "__M12_" << mlsp << ".root ";
148	>	samplename="mSUGRA";
149	>	}
150	>	if(scantype==SMS) {
151	>	copyfile << "/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "local_storage/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root ";
152	>	samplename="SMS";
153	>	}
154	>	if(scantype==GMSB) {
155	>	copyfile << "/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "local_storage/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root ";
156	>	samplename="GMSB";
157	>	}
158		stringstream newfilename;
159	<	if(scantype==mSUGRA) newfilename << "mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
160	<	if(scantype==SMS) newfilename << "SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
161	<	if(scantype==GMSB) newfilename << "GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
159	>	//    if(scantype==mSUGRA) newfilename << "local_storage/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
160	>	if(scantype==mSUGRA) newfilename << "local_storage/mSUGRA_M0_" << mglu << "__M12_" << mlsp << ".root";
161	>	if(scantype==SMS) newfilename << "local_storage/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
162	>	if(scantype==GMSB) newfilename << "local_storage/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root";
163
164		if(!Contains(((scansample.collection)[(scansample.collection).size()-1]).filename,"_"+any2string(a)+"_"+any2string(b))||(a==0&&b==0)&&!initialized_t2) {
165		dout << "The last sample is NOT the same one as the current one, possibly popping off last one and downloading as well as adding the new one." << endl;
166	<	if((scansample.collection).size()>1) {
166	>	while((scansample.collection).size()>0) {
167		scansample.RemoveLastSample();
168		}
169	+	cout << "New scanfileindex stems from scansample size: " << (scansample.collection).size() << endl;
170		scanfileindex=(scansample.collection).size();
171		//New: Loading file when necessary, not before (avoiding high memory usage and startup times)
172	<	if(scanfileindex!=0||isretry) {
173	<	delete_any_cached_scans();
174	<	dout << "Going to download the scan file with the following copy command: " << copyfile.str() << endl;
175	<	gSystem->Exec(copyfile.str().c_str());
172	>	if(1) {
173	>	//        delete_any_cached_scans();
174	>	if(!doesROOTFileExist(newfilename.str())) {
175	>	dout << "Going to download the scan file with the following copy command: " << copyfile.str() << endl;
176	>	gSystem->Exec(copyfile.str().c_str()); // download it if it hasn't been downloaded before
177	>	}
178		if(doesROOTFileExist(newfilename.str())) {
179		initialized_t2=true;
180	<	scansample.AddSample(newfilename.str(),"scansample",1,1,false,true,scanfileindex,kRed);
180	>	cout << "Going to add the following sample: " << newfilename.str() << endl;
181	>	gSystem->Exec(("ls -ltrh "+newfilename.str()).c_str());
182	>	scansample.AddSample(newfilename.str(),samplename,1,1,false,true,scanfileindex,kRed);
183		} else {
184		srmcpretries++;
185		if(srmcpretries<5) {
186		dout << "The file could not be loaded correctly - retrying!" << endl;
187		sleep(5);
188	<	load_scan_sample(a,b,scanfileindex,scantype,true);
188	>	load_scan_sample(a,b,scanfileindex,scantype,mglu,mlsp,true);
189		} else {
190		dout << "Have tried 5 times to load this sample. Giving up now and failing the program execution" << endl;
191	<	assert(0);
191	>	assert(srmcpretries<5);
192		dout << "The command " << copyfile.str() << " failed to copy the file to the local storage (saving it as " << newfilename.str() << ")" << endl;
193		}
194		}
#	Line 168 | Line 199 | void load_scan_sample(int a, int b, int
199		}
200
201		}
202	<	*/
202	>
203
204		// WATCH OUT THIS LINE NEEDS TO BE REMOVED
205		scanfileindex=0;
#	Line 205 | Line 236 | void establish_SUSY_limits(string mcjzb,
236		while ((key = (TKey*)nextkey()))
237		{
238		TObject *obj = key->ReadObj();
239	<	if(Contains((string)(obj->GetName()),"mSUGRA")) scantype=mSUGRA;
239	>	if(Contains((string)(obj->GetName()),"SUGRA")) scantype=mSUGRA;
240		if(Contains((string)(obj->GetName()),"GMSB")) scantype=GMSB;
241		}
242
#	Line 260 | Line 291 | void establish_SUSY_limits(string mcjzb,
291		TH2F *exp2mlimitmap  = new TH2F((prefix+"exp2mlimitmap"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit - 2 sigma
292
293		TH2F *exclmap  = new TH2F((prefix+"exclusionmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
263	–	TH2F *xsmap  = new TH2F((prefix+"crosssectionmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
264	–	TH2F *totxsmap  = new TH2F((prefix+"absolutecrosssectionmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
294		TH2F *flipmap  = new TH2F((prefix+"limitflipmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
295	+	TH2F *asymptoticmap  = new TH2F((prefix+"asymptoticmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
296
297
298		if(!fullerr || !mceff || !NEvents) {
#	Line 300 | Line 330 | void establish_SUSY_limits(string mcjzb,
330		vector<float> asigmas;
331		asigmas=compute_one_upper_limit(currmceff,currtoterr,ibin,mcjzb,plotfilename,true, flipped,true); // asymptotic limit first
332		float strength=asigmas[0]/XSmap->GetBinContent(GlobalBin);
333	<	if(strength>0.5&&strength<2) sigmas=compute_one_upper_limit(currmceff,currtoterr,ibin,mcjzb,plotfilename,true, flipped); // asymptotic limit first
334	<	else sigmas=asigmas;
333	>	if(strength>0.5&&strength<2) {
334	>	sigmas=compute_one_upper_limit(currmceff,currtoterr,ibin,mcjzb,plotfilename,true, flipped); // asymptotic limit first
335	>	asymptoticmap->SetBinContent(GlobalBin,0);
336	>	}
337	>	else {
338	>	sigmas=asigmas;
339	>	asymptoticmap->SetBinContent(GlobalBin,1);
340	>	}
341		exclmap->SetBinContent(GlobalBin,strength);
306	–	xsmap->SetBinContent(GlobalBin,XSmap->GetBinContent(GlobalBin));
342		}
343
309	–
310	–
311	–
344		if(sigmas[0]>0) limitmap->SetBinContent(GlobalBin,sigmas[0]); //anything else is an error code
345		if(sigmas.size()>1) {
346		explimitmap->SetBinContent(GlobalBin,sigmas[1]);
#	Line 327 | Line 359 | void establish_SUSY_limits(string mcjzb,
359		outputfile->cd();
360		limitmap->Write();
361		flipmap->Write();
362	+	asymptoticmap->Write();
363		if(doexpected) {
364		explimitmap->Write();
365		exp1plimitmap->Write();
#	Line 334 | Line 367 | void establish_SUSY_limits(string mcjzb,
367		exp2plimitmap->Write();
368		exp2mlimitmap->Write();
369		}
370	<	if(scantype==mSUGRA) {
338	<	exclmap->Write();
339	<	xsmap->Write();
340	<	totxsmap->Write();
341	<	}
370	>	if(scantype==mSUGRA) exclmap->Write();
371		outputfile->Close();
372		delete limcanvas;
373		}
#	Line 380 | Line 409 | void scan_SUSY_parameter_space(string mc
409		string prefix="SMS_";
410		// up to here, everything is set up for SMS; now we need to switch stuff around if we're dealing with an mSUGRA scan!
411		int scantype=SMS;
412	<	if(Contains((scansample.collection)[0].samplename,"mSUGRA")) scantype=mSUGRA;
412	>	if(Contains((scansample.collection)[0].samplename,"SUGRA")) scantype=mSUGRA;
413		if(Contains((scansample.collection)[0].samplename,"GMSB")) scantype=GMSB;
414
415		if(scantype==mSUGRA) {
#	Line 515 | Line 544 | write_warning(__FUNCTION__,"CURRENTLY SW
544		int a = int((PlottingSetup::ScanXzones*mlsp)/(m12trees+1));
545		int b = int((PlottingSetup::ScanYzones*mglu)/(m0trees+1));
546		int scanfileindex=PlottingSetup::ScanYzones*a+b;
547	<	load_scan_sample(a,b,scanfileindex,scantype);
547	>	load_scan_sample(a,b,scanfileindex,scantype,mglu,mlsp);
548
549		clock_t start,finish;
550		start = clock(); // starting the clock to measure how long the computation takes!
#	Line 578 | Line 607 | write_warning(__FUNCTION__,"CURRENTLY SW
607		}
608
609
610	+	if(scantype==mSUGRA) {
611	+	float absxs=0;
612	+	for(int i=0;i<12;i++) absxs+=GetXSecForPointAndChannel(mglu,mlsp,xsec,i);
613	+	XSmap->SetBinContent(GlobalBin,absxs);
614	+	}
615	+
616		if(nevents!=0&&(efficiencyonly||systematicsonly)) {
617		Value effwosigcont = MCefficiency((scansample.collection)[scanfileindex].events,result,resulterr,flipped,mcjzb,requireZ,(int)nevents,scantype,xsec,addcut.str(),-1);
618		if(result<0&&allowflipping) {
#	Line 628 | Line 663 | write_warning(__FUNCTION__,"CURRENTLY SW
663		float mcwoscefr= systematics[0][9]; // error on efficiency without signal contamination
664		float sys_pdf   = 0;
665		if(systematics[0].size()>10) sys_pdf = systematics[0][10]; // PDF
666	<
666	>	cout << "Going to store: total error is: " << toterr << " and the relative tot error, i.e. total error / mc efficiency is " << toterr/mceff << endl;
667		if(mceff!=mceff||toterr!=toterr||mceff<0 && (!systematicsonly&&!efficiencyonly)) {
668		dout << "Limits can't be calculated for this configuration (mglu="<<mglu<<" , mlsp="<<mlsp << ") as either the efficiency or its error are not positive numbers! (mceff="<<mceff<<" and toterr="<<toterr<<")"<< endl;
669		continue;

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines