- Timestamp:
- 03/08/16 15:57:17 (9 years ago)
- Location:
- trunk/Mars/mjobs
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/Mars/mjobs/MJSimulation.cc
r18333 r18450 56 56 #include <TEnv.h> 57 57 #include <TCanvas.h> 58 #include <iostream> 58 59 59 60 // Core … … 144 145 return kTRUE; 145 146 } 146 /*147 TString MJSimulation::GetOutFile(const MSequence &seq) const148 {149 return seq.IsValid() ? Form("ceres%08d.root", seq.GetSequence()) : "ceres.root";150 }151 */152 147 153 148 Bool_t MJSimulation::WriteResult(const MParList &plist, const MSequence &seq, Int_t num) … … 337 332 } 338 333 334 // Setup the header accordingly to the used operation mode 335 void MJSimulation::SetupHeaderOperationMode(MRawRunHeader &header) const 336 { 337 switch (fOperationMode) 338 { 339 case kModeData: 340 header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTData); 341 header.SetRunInfo(0, fRunNumber<0 ? 3 : fRunNumber); 342 break; 343 344 case kModePed: 345 header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTPedestal); 346 header.SetSourceInfo("Pedestal"); 347 header.SetRunInfo(0, fRunNumber<0 ? 1 : fRunNumber); 348 break; 349 350 case kModeCal: 351 header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTCalibration); 352 header.SetSourceInfo("Calibration"); 353 header.SetRunInfo(0, fRunNumber<0 ? 2 : fRunNumber); 354 break; 355 356 case kModePointRun: 357 header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTPointRun); 358 header.SetRunInfo(0, fRunNumber<0 ? 0 : fRunNumber); 359 break; 360 } 361 } 362 363 void MJSimulation::CreateBinningObjects(MParList &plist) const 364 { 365 MBinning *binse = (MBinning*) plist.FindCreateObj("MBinning","BinningEnergy"); 366 binse->SetEdges(120, 1, 1000000,"log"); 367 MBinning *binsth = (MBinning*) plist.FindCreateObj("MBinning","BinningThreshold"); 368 binsth->SetEdges( 60, 0.9, 900000, "log"); 369 MBinning *binsee = (MBinning*) plist.FindCreateObj("MBinning","BinningEnergyEst"); 370 binsee->SetEdges( 36, 0.9, 900000, "log"); 371 MBinning *binss = (MBinning*) plist.FindCreateObj("MBinning","BinningSize"); 372 binss->SetEdges( 100, 1, 10000000, "log"); 373 MBinning *binsi = (MBinning*) plist.FindCreateObj("MBinning","BinningImpact"); 374 binsi->SetEdges( 32, 0, 800); 375 MBinning *binsh = (MBinning*) plist.FindCreateObj("MBinning","BinningHeight"); 376 binsh->SetEdges( 150, 0, 50); 377 MBinning *binsaz = (MBinning*) plist.FindCreateObj("MBinning","BinningAz"); 378 binsaz->SetEdges(720, -360, 360); 379 MBinning *binszd = (MBinning*) plist.FindCreateObj("MBinning","BinningZd"); 380 binszd->SetEdges( 70, 0, 70); 381 MBinning *binsvc = (MBinning*) plist.FindCreateObj("MBinning","BinningViewCone"); 382 binsvc->SetEdges( 35, 0, 7); 383 MBinning *binsel = (MBinning*) plist.FindCreateObj("MBinning","BinningTotLength"); 384 binsel->SetEdges(150, 0, 50); 385 MBinning *binsew = (MBinning*) plist.FindCreateObj("MBinning","BinningMedLength"); 386 binsew->SetEdges(150, 0, 15); 387 plist.FindCreateObj("MBinning","BinningDistC"); 388 plist.FindCreateObj("MBinning","BinningDist"); 389 plist.FindCreateObj("MBinning","BinningTrigPos"); 390 } 391 339 392 Bool_t MJSimulation::Process(const MArgs &args, const MSequence &seq) 340 393 { 394 // -------------------------------------------------------------------------------- 395 // -------------------------------------------------------------------------------- 396 // Initialization 397 // -------------------------------------------------------------------------------- 398 // -------------------------------------------------------------------------------- 399 400 // -------------------------------------------------------------------------------- 401 // Logging output: 402 // -------------------------------------------------------------------------------- 403 // - description of the job itself 404 // - list of the to processed sequence 341 405 *fLog << inf; 342 406 fLog->Separator(GetDescriptor()); … … 357 421 return kFALSE; 358 422 } 359 360 // --------------------------------------------------------------------------------361 // Setup Parlist423 // -------------------------------------------------------------------------------- 424 // Setup MParList and MTasklist 425 // -------------------------------------------------------------------------------- 362 426 MParList plist; 363 427 plist.AddToList(this); // take care of fDisplay! … … 365 429 MTaskList tasks; 366 430 plist.AddToList(&tasks); 367 // -------------------------------------------------------------------------------- 368 431 432 // -------------------------------------------------------------------------------- 433 // -------------------------------------------------------------------------------- 434 // Parameter Container Setup 435 // -------------------------------------------------------------------------------- 436 // -------------------------------------------------------------------------------- 437 438 // -------------------------------------------------------------------------------- 439 // Setup container for the reflector, the cones and the camera geometry 440 // -------------------------------------------------------------------------------- 369 441 // FIXME: Allow empty GeomCones! 370 442 MParEnv env0("Reflector"); … … 377 449 plist.AddToList(&env1); 378 450 plist.AddToList(&env2); 379 451 // -------------------------------------------------------------------------------- 452 // Setup container for the "camera array" of the extracted number of photons 453 // -------------------------------------------------------------------------------- 454 // from ExtractorRndm 380 455 plist.FindCreateObj("MPedPhotCam", "MPedPhotFromExtractorRndm"); 381 456 // -------------------------------------------------------------------------------- 457 // Setup spline containers for: 458 // -------------------------------------------------------------------------------- 459 // - the pulse shape 460 // - the different photon acceptance splines 382 461 MParSpline shape("PulseShape"); 383 plist.AddToList(&shape);384 385 // *** FIXME *** FIXME *** FIXME ***386 plist.FindCreateObj("MMcRunHeader");387 388 MRawRunHeader header;389 header.SetValidMagicNumber();390 //header.InitFadcType(3);391 392 switch (fOperationMode)393 {394 case kModeData:395 header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTData);396 header.SetRunInfo(0, fRunNumber<0 ? 3 : fRunNumber);397 break;398 399 case kModePed:400 header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTPedestal);401 header.SetSourceInfo("Pedestal");402 header.SetRunInfo(0, fRunNumber<0 ? 1 : fRunNumber);403 break;404 405 case kModeCal:406 header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTCalibration);407 header.SetSourceInfo("Calibration");408 header.SetRunInfo(0, fRunNumber<0 ? 2 : fRunNumber);409 break;410 411 case kModePointRun:412 header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTPointRun);413 header.SetRunInfo(0, fRunNumber<0 ? 0 : fRunNumber);414 break;415 }416 417 // FIXME: Move to MSimPointingPos, MSimCalibrationSignal418 // Can we use this as input for MSimPointingPos?419 header.SetObservation("On", "MonteCarlo");420 plist.AddToList(&header);421 422 MCorsikaRead read;423 read.SetForceMode(fForceMode);424 425 if (!seq.IsValid())426 {427 for (int i=0; i<args.GetNumArguments(); i++)428 read.AddFile(args.GetArgumentStr(i));429 }430 else431 read.AddFiles(iter);432 433 MContinue precut("", "PreCut");434 precut.IsInverted();435 precut.SetAllowEmpty();436 437 MSimMMCS simmmcs;438 439 462 MParSpline splinepde("PhotonDetectionEfficiency"); 440 463 MParSpline splinemirror("MirrorReflectivity"); … … 442 465 MParSpline splinecones2("ConesTransmission"); 443 466 MParSpline splineAdditionalPhotonAcceptance("AdditionalPhotonAcceptance"); 467 plist.AddToList(&shape); 444 468 plist.AddToList(&splinepde); 445 469 plist.AddToList(&splinemirror); … … 448 472 plist.AddToList(&splineAdditionalPhotonAcceptance); 449 473 474 // -------------------------------------------------------------------------------- 475 // Setup container for the MC run header and the Raw run header 476 // -------------------------------------------------------------------------------- 477 plist.FindCreateObj("MMcRunHeader"); 478 479 MRawRunHeader header; 480 header.SetValidMagicNumber(); 481 SetupHeaderOperationMode(header); 482 // FIXME: Move to MSimPointingPos, MSimCalibrationSignal 483 // Can we use this as input for MSimPointingPos? 484 header.SetObservation("On", "MonteCarlo"); 485 plist.AddToList(&header); 486 487 // -------------------------------------------------------------------------------- 488 // Setup container for the intended pulse position and for the trigger position 489 // -------------------------------------------------------------------------------- 490 MParameterD pulpos("IntendedPulsePos"); 491 // FIXME: Set a default which could be 1/3 of the digitization window 492 // pulpos.SetVal(40); // [ns] 493 plist.AddToList(&pulpos); 494 495 MParameterD trigger("TriggerPos"); 496 trigger.SetVal(0); 497 plist.AddToList(&trigger); 498 499 500 // -------------------------------------------------------------------------------- 501 // Setup container for the fix time offset, for residual time spread and for gapd time jitter 502 // -------------------------------------------------------------------------------- 503 // Dominik and Sebastian on: fix time offsets 504 MMatrix fix_time_offsets_between_pixels_in_ns( 505 "MFixTimeOffset","titel" 506 ); 507 plist.AddToList(&fix_time_offsets_between_pixels_in_ns); 508 509 // Jens Buss on: residual time spread 510 MParameterD resTimeSpread("ResidualTimeSpread"); 511 resTimeSpread.SetVal(0.0); 512 plist.AddToList(&resTimeSpread); 513 // Jens Buss on: residual time spread 514 MParameterD gapdTimeJitter("GapdTimeJitter"); 515 gapdTimeJitter.SetVal(0.0); 516 plist.AddToList(&gapdTimeJitter); 517 518 // -------------------------------------------------------------------------------- 519 // Creation of binning objects and apply default settings 520 // -------------------------------------------------------------------------------- 521 CreateBinningObjects(plist); 522 523 // -------------------------------------------------------------------------------- 524 // -------------------------------------------------------------------------------- 525 // Tasks Setup (Reading, Absorption, Reflector) 526 // -------------------------------------------------------------------------------- 527 // -------------------------------------------------------------------------------- 528 529 // -------------------------------------------------------------------------------- 530 // Corsika read setup 531 // -------------------------------------------------------------------------------- 532 MCorsikaRead read; 533 read.SetForceMode(fForceMode); 534 535 if (!seq.IsValid()) 536 { 537 for (int i=0; i<args.GetNumArguments(); i++) 538 read.AddFile(args.GetArgumentStr(i)); 539 } 540 else 541 read.AddFiles(iter); 542 543 // -------------------------------------------------------------------------------- 544 // Precut 545 // -------------------------------------------------------------------------------- 546 MContinue precut("", "PreCut"); 547 precut.IsInverted(); 548 precut.SetAllowEmpty(); 549 550 // -------------------------------------------------------------------------------- 551 // MSimMMCS (don't know what this is doing) and calculation of the incident angle 552 // -------------------------------------------------------------------------------- 553 MSimMMCS simmmcs; 450 554 const TString sin2 = "sin(MCorsikaEvtHeader.fZd)*sin(MCorsikaRunHeader.fZdMin*TMath::DegToRad())"; 451 555 const TString cos2 = "cos(MCorsikaEvtHeader.fZd)*cos(MCorsikaRunHeader.fZdMin*TMath::DegToRad())"; … … 457 561 calcangle.SetNameParameter("IncidentAngle"); 458 562 563 // -------------------------------------------------------------------------------- 564 // Absorption tasks 565 // -------------------------------------------------------------------------------- 566 // - Atmosphere (simatm) 567 // - photon detection efficiency (absapd) 568 // - mirror reflectivity (absmir) 569 // - angular acceptance of winston cones (cones) 570 // - transmission of winston cones (cones2) 571 // - additional photon acceptance (only motivated, not measured) (additionalPhotonAcceptance) 459 572 MSimAtmosphere simatm; 460 573 MSimAbsorption absapd("SimPhotonDetectionEfficiency"); … … 470 583 additionalPhotonAcceptance.SetParName("AdditionalPhotonAcceptance"); 471 584 585 // -------------------------------------------------------------------------------- 586 // Simulating pointing position and simulating reflector 587 // -------------------------------------------------------------------------------- 472 588 MSimPointingPos pointing; 473 589 … … 477 593 // MSimStarField stars; 478 594 595 // -------------------------------------------------------------------------------- 596 // Continue tasks 597 // -------------------------------------------------------------------------------- 598 // - Processing of the current events stops, if there aren't at least two photons (see cont3) 479 599 MContinue cont1("MPhotonEvent.GetNumPhotons<1", "ContEmpty1"); 480 600 MContinue cont2("MPhotonEvent.GetNumPhotons<1", "ContEmpty2"); … … 484 604 cont3.SetAllowEmpty(); 485 605 486 // ------------------------------------------------------------------- 487 488 MBinning binse( 120, 1, 1000000, "BinningEnergy", "log"); 489 MBinning binsth( 60, 0.9, 900000, "BinningThreshold", "log"); 490 MBinning binsee( 36, 0.9, 900000, "BinningEnergyEst", "log"); 491 MBinning binss( 100, 1, 10000000, "BinningSize", "log"); 492 // MBinning binsi( 100, -500, 500, "BinningImpact"); 493 MBinning binsi( 32, 0, 800, "BinningImpact"); 494 MBinning binsh( 150, 0, 50, "BinningHeight"); 495 MBinning binsaz(720, -360, 360, "BinningAz"); 496 MBinning binszd( 70, 0, 70, "BinningZd"); 497 MBinning binsvc( 35, 0, 7, "BinningViewCone"); 498 MBinning binsel(150, 0, 50, "BinningTotLength"); 499 MBinning binsew(150, 0, 15, "BinningMedLength"); 500 MBinning binsd("BinningDistC"); 501 MBinning binsd0("BinningDist"); 502 MBinning binstr("BinningTrigPos"); 503 504 plist.AddToList(&binsee); 505 plist.AddToList(&binse); 506 plist.AddToList(&binsth); 507 plist.AddToList(&binss); 508 plist.AddToList(&binsi); 509 plist.AddToList(&binsh); 510 plist.AddToList(&binszd); 511 plist.AddToList(&binsaz); 512 plist.AddToList(&binsvc); 513 plist.AddToList(&binsel); 514 plist.AddToList(&binsew); 515 plist.AddToList(&binstr); 516 plist.AddToList(&binsd); 517 plist.AddToList(&binsd0); 518 606 // -------------------------------------------------------------------------------- 607 // -------------------------------------------------------------------------------- 608 // Tasks Setup (Cones, SiPM-Simulation, RandomPhotons, Camera, Trigger, DAQ) 609 // -------------------------------------------------------------------------------- 610 // -------------------------------------------------------------------------------- 611 612 // -------------------------------------------------------------------------------- 613 // GeomApply, SimPSF, SimGeomCam, SimCalibrationSignal 614 // -------------------------------------------------------------------------------- 615 // - MGeomApply only resizes all MGeomCam parameter container to the actual 616 // number of pixels 617 // - MSimPSF applies an additional smearing of the photons on the camera window 618 // on default it is switched of (set to 0), but can be applied by setting: 619 // MSimPSF.fSigma: <value> 620 // in the config file. 621 // Be aware, that there is also a smearing of the photons implemented in 622 // MSimReflector by setting: 623 // MReflector.SetSigmaPSF: <value> 624 // in the config file. This is at the moment done in the default config file. 625 // - MSimGeomCam identifies which pixel was hit by which photon 626 // - MSimCalibrationSignal is only used for simulated calibration runs 627 MGeomApply apply; 628 629 MSimPSF simpsf; 630 631 MSimGeomCam simgeom; 632 simgeom.SetNameGeomCam("GeomCones"); 633 634 MSimCalibrationSignal simcal; 635 simcal.SetNameGeomCam("GeomCones"); 636 637 // -------------------------------------------------------------------------------- 638 // Simulation of random photons 639 // -------------------------------------------------------------------------------- 640 // - be aware that here default values for the nsb rate and for the dark 641 // count rate are set. They will be overwritten if the values are defined 642 // in the config file 643 // - if you want to simulate a specific nsb rate you have to set the filename 644 // in the config file to an empty string and then you can specify the NSB rate: 645 // MSimRandomPhotons.FileNameNSB: 646 // MSimRandomPhotons.FrequencyNSB: 0.0 647 648 // Sky Quality Meter: 21.82M = 2.02e-4cd/m^2 649 // 1cd = 12.566 lm / 4pi sr 650 // FIXME: Simulate photons before cones and QE! 651 MSimRandomPhotons simnsb; 652 simnsb.SetFreq(5.8, 0.004); // ph/m^2/nm/sr/ns NSB, 4MHz dark count rate 653 simnsb.SetNameGeomCam("GeomCones"); 654 // FIXME: How do we handle star-light? We need the rate but we also 655 // need to process it through the mirror! 656 657 // -------------------------------------------------------------------------------- 658 // Simulation from the SiPM to the DAQ 659 // -------------------------------------------------------------------------------- 660 // - MSimAPD simulates the whole behaviour of the SiPMs: 661 // (dead time of cells, crosstalk, afterpulses) 662 // - MSimExcessNoise adds a spread on the weight of each signal in the SiPMs 663 // - MSimBundlePhotons restructured the photon list of MPhotonEvent via a look 664 // up table. At the moment the tasks is skipped, due to the fact that there 665 // is no file name specified in the config file 666 // - MSimSignalCam only fills a SignalCam container for the cherenkov photons 667 // - MSimCamera simulates the behaviour of the camera: 668 // (electronic noise, accoupling, time smearing and offset for the photons, 669 // pulses injected by all photons) 670 // - MSimTrigger simulates the behaviour of the Trigger: 671 // (Adding patch voltages, applying clipping, applying discriminator, applying 672 // time over threshold and a possible trigger logic) 673 // - MContinue conttrig stops the processing of the current event if the trigger 674 // position is negativ (and therefore not valid) 675 // - MSimReadout simulates the behaviour of the readout: 676 // (Digitization and saturation of the ADC) 677 MSimAPD simapd; 678 simapd.SetNameGeomCam("GeomCones"); 679 680 MSimExcessNoise simexcnoise; 681 MSimBundlePhotons simsum; 682 MSimSignalCam simsignal; 683 MSimCamera simcam; 684 MSimTrigger simtrig; 685 MContinue conttrig("TriggerPos.fVal<0", "ContTrigger"); 686 MSimReadout simdaq; 687 688 // -------------------------------------------------------------------------------- 689 // Standard analysis with hillas parameters for the true cherenkov photons 690 // -------------------------------------------------------------------------------- 691 // Remove isolated pixels 692 MImgCleanStd clean(0, 0); 693 //clean.SetCleanLvl0(0); // The level above which isolated pixels are kept 694 clean.SetCleanRings(0); 695 clean.SetMethod(MImgCleanStd::kAbsolute); 696 697 //clean.SetNamePedPhotCam("MPedPhotFromExtractorRndm"); 698 699 MHillasCalc hcalc; 700 hcalc.Disable(MHillasCalc::kCalcConc); 701 702 703 // -------------------------------------------------------------------------------- 704 // Setup of histogram containers and the corresponding fill tasks 705 // -------------------------------------------------------------------------------- 706 // Here is no functionality for the simulation, it is only visualization via 707 // showplot 519 708 MHn mhn1, mhn2, mhn3, mhn4; 520 709 SetupHist(mhn1); … … 565 754 plist.AddToList(&plane4); 566 755 plist.AddToList(&planeF1); 567 plist.AddToList(&planeF2); ;756 plist.AddToList(&planeF2); 568 757 569 758 //MHPSF psf; … … 589 778 //fillP.SetNameTab("PSF", "Photon distribution after cones for all mirrors"); 590 779 591 // ------------------------------------------------------------------- 592 780 781 MHCamEvent evt0a(/*10*/3, "Signal", "Average signal (absolute);;S [ph]"); 782 MHCamEvent evt0c(/*10*/3, "MaxSignal", "Maximum signal (absolute);;S [ph]"); 783 MHCamEvent evt0d(/*11*/8, "ArrTm", "Time after first photon;;T [ns]"); 784 evt0a.SetErrorSpread(kFALSE); 785 evt0c.SetCollectMax(); 786 787 MContinue cut("", "Cut"); 788 789 MFillH fillx0a(&evt0a, "MSignalCam", "FillAvgSignal"); 790 MFillH fillx0c(&evt0c, "MSignalCam", "FillMaxSignal"); 791 MFillH fillx0d(&evt0d, "MSignalCam", "FillArrTm"); 792 MFillH fillx1("MHHillas", "MHillas", "FillHillas"); 793 MFillH fillx3("MHHillasSrc", "MHillasSrc", "FillHillasSrc"); 794 MFillH fillx4("MHNewImagePar", "MNewImagePar", "FillNewImagePar"); 795 MFillH fillth("MHThreshold", "", "FillThreshold"); 796 MFillH fillca("MHCollectionArea", "", "FillTrigArea"); 797 798 fillth.SetNameTab("Threshold"); 799 fillca.SetNameTab("TrigArea"); 800 801 // -------------------------------------------------------------------------------- 802 // Formating of the outputfilepath 803 // -------------------------------------------------------------------------------- 593 804 if (!fFileOut.IsNull()) 594 805 { … … 598 809 fFileOut = fFileOut.Remove(dot); 599 810 } 600 601 // -------------------------------------------------------------------602 603 811 const char *fmt = fFileOut.IsNull() ? 604 812 Form("s/cer([0-9]+)([0-9][0-9][0-9])/%s\\/%08d.$2%%s_MonteCarlo$1.root/", Esc(fPathOut).Data(), header.GetRunNumber()) : … … 610 818 TString rule3(Form(fmt, Form("_%c", header.GetRunTypeChar()))); 611 819 820 // -------------------------------------------------------------------------------- 821 // Setup of the outputfile tasks 822 // -------------------------------------------------------------------------------- 612 823 MWriteRootFile write4a( 2, rule4, fOverwrite?"RECREATE":"NEW", "Star file"); 613 824 MWriteRootFile write4b( 2, rule4, fOverwrite?"RECREATE":"NEW", "Star file"); … … 682 893 write4b.AddContainer("MMcEvtBasic", "OriginalMC"); 683 894 684 // ------------------------------------------------------------------- 685 686 MGeomApply apply; 687 688 MSimPSF simpsf; 689 690 MSimGeomCam simgeom; 691 simgeom.SetNameGeomCam("GeomCones"); 692 693 MSimCalibrationSignal simcal; 694 simcal.SetNameGeomCam("GeomCones"); 695 696 // Sky Quality Meter: 21.82M = 2.02e-4cd/m^2 697 // 1cd = 12.566 lm / 4pi sr 698 699 // FIXME: Simulate photons before cones and QE! 700 701 MSimRandomPhotons simnsb; 702 simnsb.SetFreq(5.8, 0.004); // ph/m^2/nm/sr/ns NSB, 4MHz dark count rate 703 simnsb.SetNameGeomCam("GeomCones"); 704 705 // FIXME: How do we handle star-light? We need the rate but we also 706 // need to process it through the mirror! 707 708 MSimAPD simapd; 709 simapd.SetNameGeomCam("GeomCones"); 710 711 MSimExcessNoise simexcnoise; 712 MSimBundlePhotons simsum; 713 MSimSignalCam simsignal; 714 MSimCamera simcam; 715 MSimTrigger simtrig; 716 MSimReadout simdaq; 717 718 MContinue conttrig("TriggerPos.fVal<0", "ContTrigger"); 719 720 MParameterD pulpos("IntendedPulsePos"); 721 // FIXME: Set a default which could be 1/3 of the digitization window 722 // pulpos.SetVal(40); // [ns] 723 plist.AddToList(&pulpos); 724 725 MParameterD trigger("TriggerPos"); 726 trigger.SetVal(0); 727 plist.AddToList(&trigger); 728 729 // ------------------------------------------------------------------- 730 // Dominik and Sebastian on: fix time offsets 731 MMatrix fix_time_offsets_between_pixels_in_ns( 732 "MFixTimeOffset","titel" 733 ); 734 plist.AddToList(&fix_time_offsets_between_pixels_in_ns); 735 736 // ------------------------------------------------------------------- 737 738 // Jens Buss on: residual time spread 739 MParameterD resTimeSpread("ResidualTimeSpread"); 740 resTimeSpread.SetVal(0.0); 741 plist.AddToList(&resTimeSpread); 742 743 // ------------------------------------------------------------------- 744 745 // Jens Buss on: residual time spread 746 MParameterD gapdTimeJitter("GapdTimeJitter"); 747 gapdTimeJitter.SetVal(0.0); 748 plist.AddToList(&gapdTimeJitter); 749 750 // ------------------------------------------------------------------- 751 752 // Remove isolated pixels 753 MImgCleanStd clean(0, 0); 754 //clean.SetCleanLvl0(0); // The level above which isolated pixels are kept 755 clean.SetCleanRings(0); 756 clean.SetMethod(MImgCleanStd::kAbsolute); 757 758 //clean.SetNamePedPhotCam("MPedPhotFromExtractorRndm"); 759 760 MHillasCalc hcalc; 761 hcalc.Disable(MHillasCalc::kCalcConc); 762 763 MHCamEvent evt0a(/*10*/3, "Signal", "Average signal (absolute);;S [ph]"); 764 MHCamEvent evt0c(/*10*/3, "MaxSignal", "Maximum signal (absolute);;S [ph]"); 765 MHCamEvent evt0d(/*11*/8, "ArrTm", "Time after first photon;;T [ns]"); 766 evt0a.SetErrorSpread(kFALSE); 767 evt0c.SetCollectMax(); 768 769 MContinue cut("", "Cut"); 770 771 MFillH fillx0a(&evt0a, "MSignalCam", "FillAvgSignal"); 772 MFillH fillx0c(&evt0c, "MSignalCam", "FillMaxSignal"); 773 MFillH fillx0d(&evt0d, "MSignalCam", "FillArrTm"); 774 MFillH fillx1("MHHillas", "MHillas", "FillHillas"); 775 MFillH fillx3("MHHillasSrc", "MHillasSrc", "FillHillasSrc"); 776 MFillH fillx4("MHNewImagePar", "MNewImagePar", "FillNewImagePar"); 777 MFillH fillth("MHThreshold", "", "FillThreshold"); 778 MFillH fillca("MHCollectionArea", "", "FillTrigArea"); 779 780 fillth.SetNameTab("Threshold"); 781 fillca.SetNameTab("TrigArea"); 782 783 // ------------------------------------------------------------------- 895 896 // -------------------------------------------------------------------------------- 897 // -------------------------------------------------------------------------------- 898 // Filling of tasks in tasklist 899 // -------------------------------------------------------------------------------- 900 // -------------------------------------------------------------------------------- 784 901 785 902 if (header.IsDataRun()) 786 903 { 787 tasks.AddToList(&read); 788 tasks.AddToList(&precut); 789 tasks.AddToList(&pointing); 790 tasks.AddToList(&simmmcs); 791 if (!fPathOut.IsNull() && !HasNullOut()) 904 tasks.AddToList(&read); // Reading Corsika 905 tasks.AddToList(&precut); // Precut 906 tasks.AddToList(&pointing); // Simulating pointing 907 tasks.AddToList(&simmmcs); // Simulating MMCS 908 if (!fPathOut.IsNull() && !HasNullOut()) // Write Tasks for corsika infos 792 909 { 793 910 //tasks.AddToList(&write1b); … … 801 918 // tasks.AddToList(&stars); 802 919 if (1) 803 tasks.AddToList(&simatm); // Here because before fillh1 804 tasks.AddToList(&calcangle); 805 tasks.AddToList(&fillh1); 806 tasks.AddToList(&fillG); 920 tasks.AddToList(&simatm); // Here because before fillh1 // atmosphere absorption 921 tasks.AddToList(&calcangle); // calculation incident angle 922 tasks.AddToList(&fillh1); // fill histogram task 923 tasks.AddToList(&fillG); // fill histogram task 807 924 if (!header.IsPointRun()) 808 925 { 809 tasks.AddToList(&absapd); 810 tasks.AddToList(&absmir); 811 tasks.AddToList(&additionalPhotonAcceptance); 926 tasks.AddToList(&absapd); // photon detection efficiency of the apds 927 tasks.AddToList(&absmir); // mirror reflectivity 928 tasks.AddToList(&additionalPhotonAcceptance); // addition photon acceptance 812 929 if (0) 813 tasks.AddToList(&simatm); // FASTER? 930 tasks.AddToList(&simatm); // FASTER? // be aware this is 'commented' 814 931 } 815 tasks.AddToList(&reflect); 932 tasks.AddToList(&reflect); // Simulation of the reflector 816 933 if (!header.IsPointRun()) 817 934 { 818 tasks.AddToList(&fill0); 935 tasks.AddToList(&fill0); // fill histogram task 819 936 //tasks.AddToList(&fill1); 820 tasks.AddToList(&fill2); 821 tasks.AddToList(&fill3); 822 tasks.AddToList(&fill4); 823 tasks.AddToList(&fillF1); 937 tasks.AddToList(&fill2); // fill histogram task 938 tasks.AddToList(&fill3); // fill histogram task 939 tasks.AddToList(&fill4); // fill histogram task 940 tasks.AddToList(&fillF1); // fill histogram task 824 941 } 825 tasks.AddToList(&cones); 826 tasks.AddToList(&cones2); 942 tasks.AddToList(&cones); // angular acceptance of winston cones 943 tasks.AddToList(&cones2); // transmission of winston cones 827 944 //if (header.IsPointRun()) 828 945 // tasks.AddToList(&fillP); 829 tasks.AddToList(&cont1); 946 tasks.AddToList(&cont1); // continue if at least 2 photons 830 947 } 831 948 // ------------------------------- 832 tasks.AddToList(&apply); 949 tasks.AddToList(&apply); // apply geometry to MGeomCam containers 833 950 if (header.IsDataRun()) 834 951 { 835 tasks.AddToList(&simpsf); 836 tasks.AddToList(&simgeom); 837 tasks.AddToList(&cont2); 952 tasks.AddToList(&simpsf); // simulates additional psf (NOT used in default mode) 953 tasks.AddToList(&simgeom); // tag which pixel is hit by which photon 954 tasks.AddToList(&cont2); // continue if at least 2 photons 838 955 if (!header.IsPointRun()) 839 956 { 840 tasks.AddToList(&fillF2); 841 tasks.AddToList(&fillh2); 957 tasks.AddToList(&fillF2); // fill histogram task 958 tasks.AddToList(&fillh2); // fill histogram task 842 959 } 843 tasks.AddToList(&cont3); 960 tasks.AddToList(&cont3); // continue if at least 3 photons 844 961 } 845 962 if (fCamera) 846 963 { 847 964 if (header.IsPedestalRun() || header.IsCalibrationRun()) 848 tasks.AddToList(&simcal); 849 tasks.AddToList(&simnsb); 850 tasks.AddToList(&simapd); 851 tasks.AddToList(&simexcnoise); 852 } 853 tasks.AddToList(&simsum); 965 tasks.AddToList(&simcal); // add calibration signal for calibration runs 966 tasks.AddToList(&simnsb); // simulate nsb and dark counts 967 tasks.AddToList(&simapd); // simulate SiPM behaviour (dead time, crosstalk ...) 968 tasks.AddToList(&simexcnoise); // add excess noise 969 } 970 tasks.AddToList(&simsum); // bundle photons (NOT used in default mode) 854 971 if (fCamera) 855 972 { 856 tasks.AddToList(&simcam); 973 tasks.AddToList(&simcam); // simulate camera behaviour (creates analog signal) 857 974 if (header.IsDataRun() || fForceTrigger) 858 tasks.AddToList(&simtrig); 859 tasks.AddToList(&conttrig); 860 tasks.AddToList(&simdaq); 861 } 862 tasks.AddToList(&simsignal); // What do we do if signal<0? 975 tasks.AddToList(&simtrig); // simulate trigger 976 tasks.AddToList(&conttrig); // continue if trigger pos is valid 977 tasks.AddToList(&simdaq); // simulate data aquisition 978 } 979 tasks.AddToList(&simsignal); // What do we do if signal<0? // fill MSimSignal container 863 980 if (!fPathOut.IsNull() && !HasNullOut()) 864 981 { 865 982 //tasks.AddToList(&write1a); 866 983 if (!header.IsPedestalRun()) 867 tasks.AddToList(&write2a); 984 tasks.AddToList(&write2a); // outputtask 868 985 if (fCamera) 869 tasks.AddToList(&write3a); 986 tasks.AddToList(&write3a); // outputtask (this is the raw data writing tasks) 870 987 } 871 988 // ------------------------------- … … 874 991 // FIXME: MHCollectionArea Trigger Area! 875 992 if (header.IsDataRun()) 876 tasks.AddToList(&fillh3); 877 tasks.AddToList(&filltp); 993 tasks.AddToList(&fillh3); // fill histogram task 994 tasks.AddToList(&filltp); // fill histogram task 878 995 } 879 996 if (header.IsDataRun()) 880 997 { 881 tasks.AddToList(&fillew); 882 tasks.AddToList(&filled); 998 tasks.AddToList(&fillew); // fill histogram task 999 tasks.AddToList(&filled); // fill histogram task 883 1000 } 884 1001 if (!header.IsPedestalRun()) 885 1002 { 886 tasks.AddToList(&fillx0a); 887 tasks.AddToList(&fillx0c); 1003 tasks.AddToList(&fillx0a); // fill histogram task 1004 tasks.AddToList(&fillx0c); // fill histogram task 888 1005 } 889 1006 if (header.IsDataRun()) 890 1007 { 891 tasks.AddToList(&clean); 892 tasks.AddToList(&hcalc); 1008 tasks.AddToList(&clean); // Cleaning for standard analysis 1009 tasks.AddToList(&hcalc); // Hillas parameter calculation 893 1010 tasks.AddToList(&cut); 894 tasks.AddToList(&fillx0d); 895 tasks.AddToList(&fillx1); 1011 tasks.AddToList(&fillx0d); // fill histogram task 1012 tasks.AddToList(&fillx1); // fill histogram task 896 1013 //tasks.AddToList(&fillx2); 897 tasks.AddToList(&fillx3); 898 tasks.AddToList(&fillx4); 1014 tasks.AddToList(&fillx3); // fill histogram task 1015 tasks.AddToList(&fillx4); // fill histogram task 899 1016 if (!HasNullOut()) 900 1017 tasks.AddToList(&write4a); 901 1018 //tasks.AddToList(&fillx5); 902 1019 903 tasks.AddToList(&fillh4); 904 tasks.AddToList(&fillth); 905 tasks.AddToList(&fillca); 906 } 907 //------------------------------------------- 908 1020 tasks.AddToList(&fillh4); // fill histogram task 1021 tasks.AddToList(&fillth); // fill histogram task 1022 tasks.AddToList(&fillca); // fill histogram task 1023 } 1024 1025 1026 // -------------------------------------------------------------------------------- 1027 // -------------------------------------------------------------------------------- 1028 // Event loop and processing display 1029 // -------------------------------------------------------------------------------- 1030 // -------------------------------------------------------------------------------- 1031 1032 // -------------------------------------------------------------------------------- 1033 // Creation of event loop 1034 // -------------------------------------------------------------------------------- 909 1035 tasks.SetAccelerator(MTask::kAccDontReset|MTask::kAccDontTime); 910 1036 … … 918 1044 919 1045 // FIXME: If pedestal file given use the values from this file 920 //------------------------------------------- 1046 1047 1048 // -------------------------------------------------------------------------------- 1049 // Preparation of the graphicalk display which is shown while processing 1050 // -------------------------------------------------------------------------------- 921 1051 922 1052 MGeomCam *cam = static_cast<MGeomCam*>(env2.GetCont()); 923 1053 924 if (binstr.IsDefault()) 925 binstr.SetEdgesLin(150, -shape.GetWidth(), 1054 1055 MBinning *binsd = (MBinning*) plist.FindObject("BinningDistC"); 1056 MBinning *binsd0 = (MBinning*) plist.FindObject("BinningDist"); 1057 MBinning *binstr = (MBinning*) plist.FindObject("BinningTrigPos"); 1058 if (binstr->IsDefault()) 1059 binstr->SetEdgesLin(150, -shape.GetWidth(), 926 1060 header.GetFreqSampling()/1000.*header.GetNumSamples()+shape.GetWidth()); 927 1061 928 if (binsd .IsDefault() && cam)929 binsd .SetEdgesLin(100, 0, cam->GetMaxRadius()*cam->GetConvMm2Deg());930 931 if (binsd0 .IsDefault() && cam)932 binsd0 .SetEdgesLin(100, 0, cam->GetMaxRadius()*cam->GetConvMm2Deg());1062 if (binsd->IsDefault() && cam) 1063 binsd->SetEdgesLin(100, 0, cam->GetMaxRadius()*cam->GetConvMm2Deg()); 1064 1065 if (binsd0->IsDefault() && cam) 1066 binsd0->SetEdgesLin(100, 0, cam->GetMaxRadius()*cam->GetConvMm2Deg()); 933 1067 934 1068 … … 1045 1179 } 1046 1180 1047 //------------------------------------------- 1181 // -------------------------------------------------------------------------------- 1182 // Perform event loop 1183 // -------------------------------------------------------------------------------- 1048 1184 1049 1185 // Execute first analysis -
trunk/Mars/mjobs/MJSimulation.h
r17866 r18450 40 40 void SetupHeaderKeys(MWriteFitsFile& write, MRawRunHeader &header) const; 41 41 void SetupVetoColumns(MWriteFitsFile& write) const; 42 void SetupHeaderOperationMode(MRawRunHeader &header) const; 43 void CreateBinningObjects(MParList &plist) const; 42 44 43 45 public:
Note:
See TracChangeset
for help on using the changeset viewer.