| 1 | /* ======================================================================== *\
|
|---|
| 2 | !
|
|---|
| 3 | ! *
|
|---|
| 4 | ! * This file is part of MARS, the MAGIC Analysis and Reconstruction
|
|---|
| 5 | ! * Software. It is distributed to you in the hope that it can be a useful
|
|---|
| 6 | ! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
|
|---|
| 7 | ! * It is distributed WITHOUT ANY WARRANTY.
|
|---|
| 8 | ! *
|
|---|
| 9 | ! * Permission to use, copy, modify and distribute this software and its
|
|---|
| 10 | ! * documentation for any purpose is hereby granted without fee,
|
|---|
| 11 | ! * provided that the above copyright notice appear in all copies and
|
|---|
| 12 | ! * that both that copyright notice and this permission notice appear
|
|---|
| 13 | ! * in supporting documentation. It is provided "as is" without express
|
|---|
| 14 | ! * or implied warranty.
|
|---|
| 15 | ! *
|
|---|
| 16 | !
|
|---|
| 17 | !
|
|---|
| 18 | ! Author(s): Wolfgang Wittek 4/2002 <mailto:wittek@mppmu.mpg.de>
|
|---|
| 19 | !
|
|---|
| 20 | ! Copyright: MAGIC Software Development, 2000-2002
|
|---|
| 21 | !
|
|---|
| 22 | !
|
|---|
| 23 | \* ======================================================================== */
|
|---|
| 24 | //////////////////////////////////////////////////////////////////////////////
|
|---|
| 25 | // //
|
|---|
| 26 | // This macro calculates photon fluxes as a function of //
|
|---|
| 27 | // - energy and time //
|
|---|
| 28 | // - energy and Theta //
|
|---|
| 29 | // //
|
|---|
| 30 | // It is assumed that the data have been taken in the wobble mode. //
|
|---|
| 31 | // This means that there is only one set of data, from which both //
|
|---|
| 32 | // 'on' and 'off' data are constructed. //
|
|---|
| 33 | // //
|
|---|
| 34 | // Necessary input from MC : //
|
|---|
| 35 | // - migration matrix (E_est, E_true) as a functioin of Theta //
|
|---|
| 36 | // - effective collection areas as a function of energy and Theta //
|
|---|
| 37 | // //
|
|---|
| 38 | // //
|
|---|
| 39 | // The input from MC has to refer to the wobble mode too. //
|
|---|
| 40 | // //
|
|---|
| 41 | // The macro includes : //
|
|---|
| 42 | // - the calculation of Hillas parameters //
|
|---|
| 43 | // - the calculation of the effective on time //
|
|---|
| 44 | // - the unfolding of the distributions in the estimated energy E_est //
|
|---|
| 45 | // //
|
|---|
| 46 | // For testing purposes (as long as no real data, no energy estimator, //
|
|---|
| 47 | // no migration matrices and no effective collection areas are available) //
|
|---|
| 48 | // - event times are generated (according to a Poissonian with dead time)//
|
|---|
| 49 | // - a dummy energy estimator is provided //
|
|---|
| 50 | // - migration matrices are generated //
|
|---|
| 51 | // - a kind of effective collection areas are calculated //
|
|---|
| 52 | // //
|
|---|
| 53 | //////////////////////////////////////////////////////////////////////////////
|
|---|
| 54 | void flux()
|
|---|
| 55 | {
|
|---|
| 56 | //--------------------------------------------------------------
|
|---|
| 57 | // empty lists of parameter containers and tasks
|
|---|
| 58 | //
|
|---|
| 59 | MParList parlist;
|
|---|
| 60 | MTaskList tasklist;
|
|---|
| 61 |
|
|---|
| 62 | //--------------------------------------------------------------
|
|---|
| 63 | // Geometry information of the MAGIC camera
|
|---|
| 64 | //
|
|---|
| 65 | MGeomCamMagic geomcam;
|
|---|
| 66 |
|
|---|
| 67 | //--------------------------------------------------------------
|
|---|
| 68 | // Define the two source positions
|
|---|
| 69 | //
|
|---|
| 70 | MSrcPosCam source("Source");
|
|---|
| 71 | MSrcPosCam antisrc("AntiSource");
|
|---|
| 72 | source.SetXY(0, 0);
|
|---|
| 73 | antisrc.SetXY(+240, 0);
|
|---|
| 74 |
|
|---|
| 75 | //--------------------------------------------------------------
|
|---|
| 76 | // Setup binning for the histograms
|
|---|
| 77 | //
|
|---|
| 78 | MBinning binswidth("BinningWidth");
|
|---|
| 79 | binswidth.SetEdges(100, 0, 1); // 100 bins from 0 to 1 deg
|
|---|
| 80 |
|
|---|
| 81 | MBinning binslength("BinningLength");
|
|---|
| 82 | binslength.SetEdges(100, 0, 1); // 100 bins from 0 to 1 deg
|
|---|
| 83 |
|
|---|
| 84 | MBinning binsdist("BinningDist");
|
|---|
| 85 | binsdist.SetEdges(100, 0, 2); // 100 bins from 0 to 2 deg
|
|---|
| 86 |
|
|---|
| 87 | MBinning binsalpha("BinningAlpha");
|
|---|
| 88 | binsalpha.SetEdges(45, -90, 90);
|
|---|
| 89 |
|
|---|
| 90 | MBinning binse("BinningE");
|
|---|
| 91 | binse.SetEdgesLog(10, 10, 1e3);
|
|---|
| 92 |
|
|---|
| 93 | MBinning binstheta("BinningTheta");
|
|---|
| 94 | binstheta.SetEdges(7, -2.5, 32.5);
|
|---|
| 95 |
|
|---|
| 96 | MBinning binstime("BinningTime");
|
|---|
| 97 | binstime.SetEdges(5, 0, 100);
|
|---|
| 98 |
|
|---|
| 99 | MBinning binsdifftime("BinningTimeDiff");
|
|---|
| 100 | binsdifftime.SetEdges(50, 0, 0.1);
|
|---|
| 101 |
|
|---|
| 102 |
|
|---|
| 103 | //--------------------------------------------------------------
|
|---|
| 104 | // Fill list of parameter containers
|
|---|
| 105 | //
|
|---|
| 106 |
|
|---|
| 107 | parlist.AddToList(&tasklist);
|
|---|
| 108 | parlist.AddToList(&geomcam);
|
|---|
| 109 |
|
|---|
| 110 | parlist.AddToList(&source);
|
|---|
| 111 | parlist.AddToList(&antisrc);
|
|---|
| 112 |
|
|---|
| 113 | parlist.AddToList(&binswidth);
|
|---|
| 114 | parlist.AddToList(&binslength);
|
|---|
| 115 | parlist.AddToList(&binsdist);
|
|---|
| 116 |
|
|---|
| 117 | parlist.AddToList(&binsalpha);
|
|---|
| 118 | parlist.AddToList(&binse);
|
|---|
| 119 | parlist.AddToList(&binstheta);
|
|---|
| 120 | parlist.AddToList(&binstime);
|
|---|
| 121 | parlist.AddToList(&binsdifftime);
|
|---|
| 122 |
|
|---|
| 123 | //--------------------------------------------------------------
|
|---|
| 124 | // Setup the tasks
|
|---|
| 125 | //
|
|---|
| 126 |
|
|---|
| 127 | //.....................................
|
|---|
| 128 | // input file
|
|---|
| 129 | //
|
|---|
| 130 | // next statement commented out
|
|---|
| 131 | // MReadMarsFile reader("Events", "~/data/Gamma*.root");
|
|---|
| 132 |
|
|---|
| 133 | //MReadMarsFile reader("Events", "/hd31/rkb/Camera/mpi/Gamma*.root");
|
|---|
| 134 | //reader.AddFile("/hd31/rkb/Camera/mpi/prot_N_*.root");
|
|---|
| 135 |
|
|---|
| 136 | // next 2 statements added
|
|---|
| 137 | MReadMarsFile reader("Events");
|
|---|
| 138 | reader.AddFile("protondata/gamma_10_cn*.root");
|
|---|
| 139 |
|
|---|
| 140 | reader.EnableBranch("MMcEvt.fTheta");
|
|---|
| 141 |
|
|---|
| 142 | //.....................................
|
|---|
| 143 | // generation of event time
|
|---|
| 144 | MMcTimeGenerate rand;
|
|---|
| 145 |
|
|---|
| 146 | //.....................................
|
|---|
| 147 | // collection area
|
|---|
| 148 | MMcCollectionAreaCalc acalc;
|
|---|
| 149 |
|
|---|
| 150 | MMcPedestalCopy pcopy;
|
|---|
| 151 | MMcPedestalNSBAdd pnsb;
|
|---|
| 152 | MCerPhotCalc ncalc;
|
|---|
| 153 | MImgCleanStd clean;
|
|---|
| 154 | MBlindPixelCalc blind;
|
|---|
| 155 |
|
|---|
| 156 | //.....................................
|
|---|
| 157 | // source independent image parameters
|
|---|
| 158 | MHillasCalc hcalc;
|
|---|
| 159 |
|
|---|
| 160 | //.....................................
|
|---|
| 161 | // source dependent image parameters
|
|---|
| 162 | MHillasSrcCalc hsrc1("Source", "HillasSrc");
|
|---|
| 163 | MHillasSrcCalc hsrc2("AntiSource", "HillasAntiSrc");
|
|---|
| 164 |
|
|---|
| 165 | //.....................................
|
|---|
| 166 | // energy estimation
|
|---|
| 167 | MEnergyEstimate estim;
|
|---|
| 168 |
|
|---|
| 169 | //.....................................
|
|---|
| 170 | // migration matrix (MC)
|
|---|
| 171 | MFillH eestetrue("MHMcEnergyMigration", "HillasSrc");
|
|---|
| 172 |
|
|---|
| 173 |
|
|---|
| 174 | //.....................................
|
|---|
| 175 | // plots of width and length
|
|---|
| 176 | MFillH hfill1h("MHHillas", "MHillas");
|
|---|
| 177 |
|
|---|
| 178 | //.....................................
|
|---|
| 179 | // star map
|
|---|
| 180 | MFillH hfill1m("MHStarMap", "MHillas");
|
|---|
| 181 |
|
|---|
| 182 | //.....................................
|
|---|
| 183 | // plots of alpha and dist
|
|---|
| 184 | MFillH hfill2s("HSource [MHHillasSrc]", "HillasSrc");
|
|---|
| 185 | MFillH hfill2a("HAntiSource [MHHillasSrc]", "HillasAntiSrc");
|
|---|
| 186 |
|
|---|
| 187 | //.....................................
|
|---|
| 188 | // filters
|
|---|
| 189 | MFTriggerLvl1 lvl1;
|
|---|
| 190 |
|
|---|
| 191 | const Float_t alpha0 = 15; // [deg]
|
|---|
| 192 | MFAlpha fsrc ("HillasSrc", '>', alpha0);
|
|---|
| 193 | MFAlpha fasrc("HillasAntiSrc", '>', alpha0);
|
|---|
| 194 |
|
|---|
| 195 | //.....................................
|
|---|
| 196 | // 1-D profile plots (<Theta>, time)
|
|---|
| 197 | // (<Theta>, Theta)
|
|---|
| 198 | MFillH fillthetabartime ("ThetabarTime [MHThetabarTime]", "MMcEvt");
|
|---|
| 199 | // fillthetabartime.SetFilter(&lvl1);
|
|---|
| 200 | MFillH fillthetabartheta("ThetabarTheta [MHThetabarTheta]", "MMcEvt");
|
|---|
| 201 | // fillthetabartheta.SetFilter(&lvl1);
|
|---|
| 202 |
|
|---|
| 203 | //.....................................
|
|---|
| 204 | // 2-D plots (Delta(time), time)
|
|---|
| 205 | // (Delta(time), Theta)
|
|---|
| 206 | MFillH filldtimetime ("EffOnTime [MHTimeDiffTime]", "MMcEvt");
|
|---|
| 207 | // filldtime.SetFilter(&lvl1);
|
|---|
| 208 | MFillH filldtimetheta("EffOnTheta [MHTimeDiffTheta]", "MMcEvt");
|
|---|
| 209 | // fillontheta.SetFilter(&lvl1);
|
|---|
| 210 |
|
|---|
| 211 | //.....................................
|
|---|
| 212 | // 3-D plots (alpha, E_est, time)
|
|---|
| 213 | MFillH fillsptime ("FluxSrcTime [MHAlphaEnergyTime]", "HillasSrc");
|
|---|
| 214 | fillsptime.SetFilter(&fasrc);
|
|---|
| 215 | MFillH fillasptime ("FluxASrcTime [MHAlphaEnergyTime]", "HillasAntiSrc");
|
|---|
| 216 | fillasptime.SetFilter(&fsrc);
|
|---|
| 217 |
|
|---|
| 218 | //.....................................
|
|---|
| 219 | // 3-D plots (alpha, E_est, Theta)
|
|---|
| 220 | MFillH fillsptheta ("FluxSrcTheta [MHAlphaEnergyTheta]", "HillasSrc");
|
|---|
| 221 | fillsptheta.SetFilter(&fasrc);
|
|---|
| 222 | MFillH fillasptheta("FluxASrcTheta [MHAlphaEnergyTheta]", "HillasAntiSrc");
|
|---|
| 223 | fillasptheta.SetFilter(&fsrc);
|
|---|
| 224 |
|
|---|
| 225 | //.....................................
|
|---|
| 226 | // MFillH fethetaall("AllTheta [MHEnergyTheta]", "MMcEvt");
|
|---|
| 227 | // MFillH fethetasel("SelTheta [MHEnergyTheta]", "MMcEvt");
|
|---|
| 228 | // fethetasel.SetFilter(&lvl1);
|
|---|
| 229 |
|
|---|
| 230 | // MFillH fetimeall ("AllTime [MHEnergyTime]", "MMcEvt");
|
|---|
| 231 | // MFillH fetimesel ("SelTime [MHEnergyTime]", "MMcEvt");
|
|---|
| 232 | // fetimesel.SetFilter(&lvl1);
|
|---|
| 233 |
|
|---|
| 234 |
|
|---|
| 235 | //---------------------------------------------------------------------
|
|---|
| 236 | // Setup the task list
|
|---|
| 237 | //
|
|---|
| 238 | tasklist.AddToList(&reader);
|
|---|
| 239 | tasklist.AddToList(&rand);
|
|---|
| 240 |
|
|---|
| 241 | tasklist.AddToList(&acalc);
|
|---|
| 242 |
|
|---|
| 243 | tasklist.AddToList(&pcopy);
|
|---|
| 244 | tasklist.AddToList(&pnsb);
|
|---|
| 245 | tasklist.AddToList(&ncalc);
|
|---|
| 246 | tasklist.AddToList(&clean);
|
|---|
| 247 | tasklist.AddToList(&blind);
|
|---|
| 248 | tasklist.AddToList(&hcalc);
|
|---|
| 249 | tasklist.AddToList(&hsrc1);
|
|---|
| 250 | tasklist.AddToList(&hsrc2);
|
|---|
| 251 | tasklist.AddToList(&estim);
|
|---|
| 252 | tasklist.AddToList(&eestetrue);
|
|---|
| 253 |
|
|---|
| 254 | tasklist.AddToList(&hfill1h);
|
|---|
| 255 | tasklist.AddToList(&hfill1m);
|
|---|
| 256 | tasklist.AddToList(&hfill2s);
|
|---|
| 257 | tasklist.AddToList(&hfill2a);
|
|---|
| 258 |
|
|---|
| 259 | tasklist.AddToList(&fillthetabartime);
|
|---|
| 260 | tasklist.AddToList(&fillthetabartheta);
|
|---|
| 261 |
|
|---|
| 262 | tasklist.AddToList(&filldtimetime);
|
|---|
| 263 | tasklist.AddToList(&filldtimetheta);
|
|---|
| 264 |
|
|---|
| 265 | // tasklist.AddToList(&lvl1);
|
|---|
| 266 |
|
|---|
| 267 | // tasklist.AddToList(&fethetaall);
|
|---|
| 268 | // tasklist.AddToList(&fethetasel);
|
|---|
| 269 | // tasklist.AddToList(&fetimeall);
|
|---|
| 270 | // tasklist.AddToList(&fetimesel);
|
|---|
| 271 |
|
|---|
| 272 | tasklist.AddToList(&fsrc);
|
|---|
| 273 | tasklist.AddToList(&fasrc);
|
|---|
| 274 |
|
|---|
| 275 | tasklist.AddToList(&fillsptime);
|
|---|
| 276 | tasklist.AddToList(&fillasptime);
|
|---|
| 277 |
|
|---|
| 278 | tasklist.AddToList(&fillsptheta);
|
|---|
| 279 | tasklist.AddToList(&fillasptheta);
|
|---|
| 280 |
|
|---|
| 281 |
|
|---|
| 282 | //----------------------------------------------------------------------
|
|---|
| 283 | // Event loop
|
|---|
| 284 | //
|
|---|
| 285 | MEvtLoop magic;
|
|---|
| 286 | magic.SetParList(&parlist);
|
|---|
| 287 |
|
|---|
| 288 | //
|
|---|
| 289 | // loop over all events
|
|---|
| 290 | //
|
|---|
| 291 | if (!magic.Eventloop())
|
|---|
| 292 | return;
|
|---|
| 293 | //----------------------------------------------------------------------
|
|---|
| 294 |
|
|---|
| 295 |
|
|---|
| 296 | tasklist.PrintStatistics(10);
|
|---|
| 297 |
|
|---|
| 298 | /*
|
|---|
| 299 | parlist.FindObject("HSource")->DrawClone();;
|
|---|
| 300 | parlist.FindObject("MHHillas")->DrawClone();
|
|---|
| 301 | */
|
|---|
| 302 | /*
|
|---|
| 303 | parlist.FindObject("MHStarMap")->DrawClone();
|
|---|
| 304 | */
|
|---|
| 305 | /*
|
|---|
| 306 | parlist.FindObject("HAntiSource")->DrawClone();
|
|---|
| 307 | */
|
|---|
| 308 | /*
|
|---|
| 309 | parlist.FindObject("MHMcCollectionArea")->DrawClone();
|
|---|
| 310 | */
|
|---|
| 311 |
|
|---|
| 312 | //----------------------------------------------------------------------
|
|---|
| 313 | // average Theta versus time
|
|---|
| 314 | // and versus Theta
|
|---|
| 315 | //
|
|---|
| 316 | MHThetabarTime *bartime = (MHThetabarTime*)parlist.FindObject("ThetabarTime");
|
|---|
| 317 | MHThetabarTheta *bartheta = (MHThetabarTheta*)parlist.FindObject("ThetabarTheta");
|
|---|
| 318 |
|
|---|
| 319 | /*
|
|---|
| 320 | bartime->DrawClone();
|
|---|
| 321 | bartheta->DrawClone();
|
|---|
| 322 | */
|
|---|
| 323 |
|
|---|
| 324 |
|
|---|
| 325 | //----------------------------------------------------------------------
|
|---|
| 326 | // this is something like the collection area
|
|---|
| 327 | // as a function of time
|
|---|
| 328 | // and as a function of Theta
|
|---|
| 329 | /*
|
|---|
| 330 | MHEnergyTime &alltime = *(MHEnergyTime*)parlist.FindObject("AllTime");
|
|---|
| 331 | MHEnergyTime &seltime = *(MHEnergyTime*)parlist.FindObject("SelTime");
|
|---|
| 332 | MHEnergyTime collareatime;
|
|---|
| 333 | collareatime.Divide(&seltime, &alltime);
|
|---|
| 334 |
|
|---|
| 335 | MHEnergyTheta &alltheta = *(MHEnergyTheta*)parlist.FindObject("AllTheta");
|
|---|
| 336 | MHEnergyTheta &seltheta = *(MHEnergyTheta*)parlist.FindObject("SelTheta");
|
|---|
| 337 | MHEnergyTheta collareatheta;
|
|---|
| 338 | collareatheta.Divide(&seltheta, &alltheta);
|
|---|
| 339 | */
|
|---|
| 340 |
|
|---|
| 341 | //----------------------------------------------------------------------
|
|---|
| 342 | // Effective on time
|
|---|
| 343 |
|
|---|
| 344 | //....................................
|
|---|
| 345 | // get plots of time differences for different intervals in time
|
|---|
| 346 | // and plots of time differences for different intervals in Theta
|
|---|
| 347 |
|
|---|
| 348 | MHTimeDiffTime *dtimetime = (MHTimeDiffTime*)parlist.FindObject("EffOnTime");
|
|---|
| 349 | MHTimeDiffTheta *dtimetheta = (MHTimeDiffTheta*)parlist.FindObject("EffOnTheta");
|
|---|
| 350 | /*
|
|---|
| 351 | dtimetime->DrawClone();
|
|---|
| 352 | dtimetheta->DrawClone();
|
|---|
| 353 | */
|
|---|
| 354 |
|
|---|
| 355 |
|
|---|
| 356 | //....................................
|
|---|
| 357 | // calculate effective on time for different intervals in time
|
|---|
| 358 | // and for different intervals in Theta
|
|---|
| 359 | MHEffOnTimeTime effontime;
|
|---|
| 360 | MHEffOnTimeTheta effontheta;
|
|---|
| 361 | effontime.SetupFill(&parlist);
|
|---|
| 362 | effontheta.SetupFill(&parlist);
|
|---|
| 363 |
|
|---|
| 364 | effontime.Calc(dtimetime->GetHist());
|
|---|
| 365 | effontheta.Calc(dtimetheta->GetHist());
|
|---|
| 366 |
|
|---|
| 367 | // plot effective on time versus time
|
|---|
| 368 | // and versus Theta
|
|---|
| 369 |
|
|---|
| 370 | effontime.DrawClone();
|
|---|
| 371 | effontheta.DrawClone();
|
|---|
| 372 |
|
|---|
| 373 |
|
|---|
| 374 |
|
|---|
| 375 | //----------------------------------------------------------------------
|
|---|
| 376 | // 3D-plots of alpha, E_est and time
|
|---|
| 377 | // and of alpha, E_est and Theta
|
|---|
| 378 | // with alpha calculated relative to the source position
|
|---|
| 379 | // and relative to the anti-source position
|
|---|
| 380 |
|
|---|
| 381 | MHAlphaEnergyTime *evtsptime = (MHAlphaEnergyTime*)parlist.FindObject("FluxSrcTime");
|
|---|
| 382 | MHAlphaEnergyTime *evtasptime = (MHAlphaEnergyTime*)parlist.FindObject("FluxASrcTime");
|
|---|
| 383 | MHAlphaEnergyTheta *evtsptheta = (MHAlphaEnergyTheta*)parlist.FindObject("FluxSrcTheta");
|
|---|
| 384 | MHAlphaEnergyTheta *evtasptheta = (MHAlphaEnergyTheta*)parlist.FindObject("FluxASrcTheta");
|
|---|
| 385 |
|
|---|
| 386 | /* this test shows that the addresses
|
|---|
| 387 | 'evtsptime' and '&evtspobj' are identical
|
|---|
| 388 | MHAlphaEnergyTime &evtspobj = *(MHAlphaEnergyTime*)parlist.FindObject("FluxSrcTime");
|
|---|
| 389 | cout << "evtsptime = " << evtsptime << endl;
|
|---|
| 390 | cout << "&evtspobj = " << &evtspobj << endl;
|
|---|
| 391 | */
|
|---|
| 392 |
|
|---|
| 393 | /*
|
|---|
| 394 | evtsptime->SetTitle("3D-plot Al-En-time, alpha\_{asp} .gt. alpha\_0");
|
|---|
| 395 | evtsptime->DrawClone();
|
|---|
| 396 | evtasptime->SetTitle("3D-plot Al-En-time, alpha\_{sp} .gt. alpha\_0");
|
|---|
| 397 | evtasptime->DrawClone();
|
|---|
| 398 |
|
|---|
| 399 | evtsptheta->SetTitle("3D-plot Al-En-Theta, alpha\_{asp} .gt. alpha\_0");
|
|---|
| 400 | evtsptheta->DrawClone();
|
|---|
| 401 | evtasptheta->SetTitle("3D-plot Al-En-Theta, alpha\_{sp} .gt. alpha\_0");
|
|---|
| 402 | evtasptheta->DrawClone();
|
|---|
| 403 | */
|
|---|
| 404 |
|
|---|
| 405 |
|
|---|
| 406 | //----------------------------------------------------------------------
|
|---|
| 407 | // Difference between source and antisource (= 'gamma' sample)
|
|---|
| 408 | // for 3D-plots of alpha, E_est and time
|
|---|
| 409 | // and of alpha, E_est and Theta
|
|---|
| 410 |
|
|---|
| 411 | MHAlphaEnergyTime gammatime("Al-En-time","3D-plot Al-En-time");
|
|---|
| 412 | MHAlphaEnergyTheta gammatheta("Al-En-Theta","3D-plot Al-En-Theta");
|
|---|
| 413 |
|
|---|
| 414 | gammatime.SetupFill(&parlist);
|
|---|
| 415 | gammatime.Subtract(evtsptime, evtasptime);
|
|---|
| 416 | gammatheta.SetupFill(&parlist);
|
|---|
| 417 | gammatheta.Subtract(evtsptheta, evtasptheta);
|
|---|
| 418 |
|
|---|
| 419 |
|
|---|
| 420 | /*
|
|---|
| 421 | gammatime.SetTitle("3D-plot Al-En-time for \'gamma\' sample");
|
|---|
| 422 | gammatime.DrawClone();
|
|---|
| 423 |
|
|---|
| 424 | gammatheta.SetTitle("3D-plot Al-En-Theta for \'gamma\' sample");
|
|---|
| 425 | gammatheta.DrawClone();
|
|---|
| 426 | */
|
|---|
| 427 |
|
|---|
| 428 | //----------------------------------------------------------------------
|
|---|
| 429 | // get number of 'gammas' in the alpha range (lo, up) as a function of
|
|---|
| 430 | // E_est and time
|
|---|
| 431 | // or E_est and Theta
|
|---|
| 432 |
|
|---|
| 433 | Axis_t lo = -10;
|
|---|
| 434 | Axis_t up = 10;
|
|---|
| 435 |
|
|---|
| 436 | // TH2D &evttime = *gammatime.GetAlphaProjection (lo, up);
|
|---|
| 437 | // TH2D &evttheta = *gammatheta.GetAlphaProjection(lo, up);
|
|---|
| 438 |
|
|---|
| 439 | TH2D *evttime = gammatime.DrawAlphaProjection (lo, up);
|
|---|
| 440 | TH2D *evttheta = gammatheta.DrawAlphaProjection(lo, up);
|
|---|
| 441 |
|
|---|
| 442 |
|
|---|
| 443 | //----------------------------------------------------------------------
|
|---|
| 444 | // plot migration matrix: E-true E-est Theta
|
|---|
| 445 |
|
|---|
| 446 | parlist.FindObject("MHMcEnergyMigration")->DrawClone();
|
|---|
| 447 |
|
|---|
| 448 | return;
|
|---|
| 449 |
|
|---|
| 450 | //----------------------------------------------------------------------
|
|---|
| 451 | for (int i=1; i<=binstime.GetNumBins(); i++)
|
|---|
| 452 | {
|
|---|
| 453 | if (effontime.GetHist()->GetBinContent(i)==0)
|
|---|
| 454 | continue;
|
|---|
| 455 |
|
|---|
| 456 | TH1D &hist = *evttime.ProjectionX("Number of Gammas vs. EnergyEst for a fixed time", i, i);
|
|---|
| 457 |
|
|---|
| 458 | /* UNFOLDING */
|
|---|
| 459 |
|
|---|
| 460 | //hist->Divide(collareatime);
|
|---|
| 461 | hist.Scale(1./effontime.GetHist()->GetBinContent(i));
|
|---|
| 462 |
|
|---|
| 463 | for (int j=1; j<=binse.GetNumBins(); j++)
|
|---|
| 464 | hist.SetBinContent(j, hist.GetBinContent(j)/hist.GetBinWidth(j));
|
|---|
| 465 |
|
|---|
| 466 | hist.SetName("Flux");
|
|---|
| 467 | hist.SetTitle("Flux[Gammas/s/m^2/GeV] vs. EnergyTrue for a fixed Time");
|
|---|
| 468 |
|
|---|
| 469 | char n[100];
|
|---|
| 470 | sprintf(n, "Canv%d", j);
|
|---|
| 471 | c= new TCanvas(n, "Title");
|
|---|
| 472 | /*
|
|---|
| 473 | hist.DrawCopy();
|
|---|
| 474 | */
|
|---|
| 475 | }
|
|---|
| 476 |
|
|---|
| 477 | delete &evttime;
|
|---|
| 478 | delete &evttheta;
|
|---|
| 479 |
|
|---|
| 480 | return;
|
|---|
| 481 |
|
|---|
| 482 | // ------------------------------------------
|
|---|
| 483 |
|
|---|
| 484 | MHMcCollectionArea carea;
|
|---|
| 485 | TH1D *collareatime = carea.GetHist(); // FIXME!
|
|---|
| 486 | TH1D *collareatheta = carea.GetHist(); // FIXME!
|
|---|
| 487 | }
|
|---|
| 488 | //===========================================================================
|
|---|
| 489 |
|
|---|
| 490 |
|
|---|
