Changeset 4234 for trunk/MagicSoft/Mars/mtemp/mifae
- Timestamp:
- 05/28/04 11:56:06 (20 years ago)
- Location:
- trunk/MagicSoft/Mars/mtemp/mifae/programs
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Mars/mtemp/mifae/programs/optimizeCuts.cc
r4098 r4234 36 36 Float_t lowdistcut=0; 37 37 Float_t uppdistcut=10; 38 Float_t sizecut=1800; 38 Float_t lowsizecut=1800; 39 Float_t uppsizecut=1800; 40 Int_t nsizebins=1; 39 41 UInt_t onRate=1; 40 42 Float_t lowlgcut=0.01; … … 130 132 char evecut[256]; 131 133 sprintf(evecut,"event%%%d==0",onRate); 132 sprintf(gencut," size>%f && dist>%f && dist<%f",sizecut,lowdistcut,uppdistcut);134 sprintf(gencut,"dist>%f && dist<%f",lowdistcut,uppdistcut); 133 135 cout << "General cut " << gencut << " (" << evecut << " for on-data)" << endl; 134 136 135 // loop on widht and length cuts 136 Float_t length=lowlgcut; 137 while(length<=upplgcut) 138 { 139 Float_t width=lowwdcut; 140 while(width<=uppwdcut) 141 { 142 TH1F* onhisto = new TH1F("OnHist" ,"Alpha Plot",nbins,0,90); 143 TH1F* offhisto = new TH1F("OffHist","Alpha Plot",nbins,0,90); 144 145 // define the width/length cut 146 char varcut[256]; 147 sprintf(varcut,"length<%f && width<%f",length,width); 148 cout << "Cutting " << varcut << endl; 149 150 // define the on/off data cut 151 char offcut[1024]; 152 sprintf(offcut,"%s && %s",gencut,varcut); 153 char oncut[1024]; 154 sprintf(oncut,"%s && %s",evecut,offcut); 155 156 // Fill the histos 157 ton->Draw("alpha>>OnHist",oncut); 158 toff->Draw("alpha>>OffHist",offcut); 159 160 // Normalize histos 161 const Int_t inibin = (Int_t)(frontiere/90.*nbins+1); 162 Float_t level=0; 163 Float_t leveloff=0; 164 for(Int_t ibin = inibin; ibin<=nbins;ibin++) 137 Bool_t isDiff = kFALSE; 138 if(nsizebins<0) 139 { 140 nsizebins*=-1; 141 isDiff=kTRUE; 142 } 143 else 144 nsizebins+=1; 145 146 // loop on size, width and length cuts 147 for(Int_t isb=0;isb<nsizebins;isb++) // loop on size 148 { 149 cout << isb << endl; 150 Float_t minsize = isDiff ? 151 TMath::Power(10,TMath::Log10(lowsizecut)+isb*(TMath::Log10(uppsizecut)-TMath::Log10(lowsizecut))/nsizebins) : 152 TMath::Power(10,TMath::Log10(lowsizecut)+isb*(TMath::Log10(uppsizecut)-TMath::Log10(lowsizecut))/(nsizebins-1)); 153 Float_t maxsize = isDiff ? TMath::Power(10,TMath::Log10(lowsizecut)+(isb+1)*(TMath::Log10(uppsizecut)-TMath::Log10(lowsizecut))/nsizebins) : 9999999999.; 154 155 Float_t length=lowlgcut; 156 while(length<=upplgcut) // loop on length 157 { 158 Float_t width=lowwdcut; 159 while(width<=uppwdcut) // loop on width 165 160 { 166 level+=onhisto->GetBinContent(ibin); 167 leveloff+=offhisto->GetBinContent(ibin); 161 TH1F* onhisto = new TH1F("OnHist" ,"Alpha Plot",nbins,0,90); 162 TH1F* offhisto = new TH1F("OffHist","Alpha Plot",nbins,0,90); 163 164 // define the width/length cut 165 char varcut[256]; 166 sprintf(varcut,"size>%f && size<%f && length<%f && width<%f",minsize,maxsize,length,width); 167 cout << "Cutting " << varcut << endl; 168 169 // define the on/off data cut 170 char offcut[1024]; 171 sprintf(offcut,"%s && %s",gencut,varcut); 172 char oncut[1024]; 173 sprintf(oncut,"%s && %s",evecut,offcut); 174 175 // Fill the histos 176 ton->Draw("alpha>>OnHist",oncut); 177 toff->Draw("alpha>>OffHist",offcut); 178 179 // Normalize histos 180 const Int_t inibin = (Int_t)(frontiere/90.*nbins+1); 181 Float_t level=0; 182 Float_t leveloff=0; 183 for(Int_t ibin = inibin; ibin<=nbins;ibin++) 184 { 185 level+=onhisto->GetBinContent(ibin); 186 leveloff+=offhisto->GetBinContent(ibin); 187 } 188 offhisto->Sumw2(); // needed to compute correct errors after normalization 189 const Float_t norm = leveloff ? level/leveloff: 1; 190 cout << "Normalizing by factor " << norm <<endl; 191 offhisto->Scale(norm); 192 193 // compute significance/excess 194 Float_t sig=0,bg=0,esig=0,ebg=0; 195 Float_t significance=0; 196 const Int_t signbins = inibin-1; 197 for(Int_t ibin = 1; ibin<=signbins;ibin++) 198 { 199 sig += onhisto->GetBinContent(ibin); 200 esig += onhisto->GetBinError(ibin)*onhisto->GetBinError(ibin); 201 bg += offhisto->GetBinContent(ibin); 202 ebg += offhisto->GetBinError(ibin)*offhisto->GetBinError(ibin); 203 } 204 Float_t error= TMath::Sqrt(esig+ebg); 205 significance = error ? (sig-bg)/error : 0; 206 207 cout << "Excess: " << sig-bg << endl; 208 cout << "N bkg: " << bg << endl; 209 cout << "Significance: " << significance << endl; 210 211 // save the result in file 212 fout << minsize << '\t' 213 << maxsize << '\t' 214 << width << '\t' 215 << length << '\t' 216 << sig-bg << '\t' 217 << significance << '\t' << endl; 218 219 delete onhisto; 220 delete offhisto; 221 222 width+=wdstep; 168 223 } 169 offhisto->Sumw2(); // needed to compute correct errors after normalization 170 const Float_t norm = leveloff ? level/leveloff: 1; 171 cout << "Normalizing by factor " << norm <<endl; 172 offhisto->Scale(norm); 173 174 // compute significance/excess 175 Float_t sig=0,bg=0,esig=0,ebg=0; 176 Float_t significance=0; 177 const Int_t signbins = inibin-1; 178 for(Int_t ibin = 1; ibin<=signbins;ibin++) 179 { 180 sig += onhisto->GetBinContent(ibin); 181 esig += onhisto->GetBinError(ibin)*onhisto->GetBinError(ibin); 182 bg += offhisto->GetBinContent(ibin); 183 ebg += offhisto->GetBinError(ibin)*offhisto->GetBinError(ibin); 184 } 185 Float_t error= TMath::Sqrt(esig+ebg); 186 significance = error ? (sig-bg)/error : 0; 187 188 cout << "Excess: " << sig-bg << endl; 189 cout << "N bkg: " << bg << endl; 190 cout << "Significance: " << significance << endl; 191 192 // save the result in file 193 fout << width << '\t' 194 << length << '\t' 195 << sig-bg << '\t' 196 << significance << '\t' << endl; 197 198 delete onhisto; 199 delete offhisto; 200 201 width+=wdstep; 202 } 203 length+=lgstep; 204 } 205 224 length+=lgstep; 225 } 226 } 206 227 fout.close(); 207 228 } … … 260 281 } 261 282 262 // size cut 263 if(strcmp(word.Data(),"SIZECUT")==0) 264 ifun >> sizecut; 283 // size cuts 284 if(strcmp(word.Data(),"SIZECUTS")==0) 285 { 286 ifun >> lowsizecut; 287 ifun >> uppsizecut; 288 ifun >> nsizebins; 289 } 265 290 266 291 // dist cuts … … 317 342 cout << "Off-data input file name(s): " << offFile << endl; 318 343 cout << "Output file name: " << outputFile << endl; 319 cout << "Size cut: " << sizecut << endl;320 344 cout << "Dist cuts (deg): [" << lowdistcut<< ","<<uppdistcut<< "]" << endl; 345 cout << "Scanning size range: [" << lowsizecut << "," << uppsizecut << "]"<< " in " << nsizebins << " steps" << endl; 321 346 cout << "Scanning length range (deg): [" << lowlgcut<< ","<<upplgcut<< "], with step (deg): "<< lgstep << endl; 322 347 cout << "Scanning width range (deg): [" << lowwdcut<< ","<<uppwdcut<< "], with step (deg): "<< wdstep << endl; -
trunk/MagicSoft/Mars/mtemp/mifae/programs/optimizecuts.datacard
r4117 r4234 1 ///////////////////////////////////////////////////////// 2 // Data card template for optimizeCuts executable // 3 // // 4 // *** Compulsory cards are: // 5 // ONFILES, OFFFILES, OUTFILE // 6 // // 7 // the rest are optional // 8 // *** Do NEVER add a datacard with no value after it // 9 ///////////////////////////////////////////////////////// 1 10 2 11 // Maximun number of on and off events to be processed) … … 4 13 5 14 // On data acceptance rate (e.g 4 for taking 1/4 of whole data sample) 6 ONRATE 215 ONRATE 3 7 16 8 17 // Input file name pattern (On data) 9 ONFILES /mnt/users/jrico/magic/mars/Mars_Standard02/mtemp/mifae/hillas/ crab20040215OnRotateCalA-D.root18 ONFILES /mnt/users/jrico/magic/mars/Mars_Standard02/mtemp/mifae/hillas/mrk20040215OnRotateNoCalB.root 10 19 11 20 // Input file name pattern (Off data) 12 OFFFILES /mnt/users/jrico/magic/mars/Mars_Standard02/mtemp/mifae/hillas/ crab20040215OffRotateCalA-H.root21 OFFFILES /mnt/users/jrico/magic/mars/Mars_Standard02/mtemp/mifae/hillas/mrk20040215OffRotateA-C.root 13 22 14 23 // output file name 15 OUTFILE ./ optimizeCrab_cal_fine.out24 OUTFILE ./prueba.out 16 25 17 // Preliminar cuts (size in size units, distance in deg)18 SIZECUT 200019 DISTCUTS 0.2 1.120 26 21 // Length initial, final and step values22 LENGTHCUTS 0.15 0.35 0.005 27 // Preliminary distance lower and upper cuts (in deg) 28 DISTCUTS 0.3 1.2 23 29 24 // Width initial, final and step values 25 WIDTHCUTS 0.05 0.15 0.005 30 // Size initial, final and number of steps 31 // it will take logarithmic constant intervals 32 // positive number of steps means lower cut (integral) 33 // negative number of steps means lower and upper cuts (diferential) 34 SIZECUTS 500 10000 3 35 36 // Length initial, final and step values (in deg) 37 LENGTHCUTS 0.15 0.35 0.1 38 39 // Width initial, final and step values (in deg) 40 WIDTHCUTS 0.05 0.15 0.1 26 41 27 42
Note:
See TracChangeset
for help on using the changeset viewer.