Notes and Further References

1. Start by looking at the state UserAnalysis has by looking at its data members which will be found in the private: section of UserAnalysis.h:- private: TH1F fQtotHisto; // Histogram of total charge TH1F fNhitHisto; // Histogram of total number of hit strips TH1F fEvtQHisto; // Histogram of charges for a single data record TH1F fEvtTHisto; // Histogram of times for a single data record TNtuple fNtuple; // Ntuple of event summary data int fNevents; // Number of events seen double fSumQ; // Accumulated sum of charge double fSumQsqr; // Accumulated sum of charge^2 double fQminCut; // A lower bound on digit charge double fQtotMaxCut; // A cut value to make pass/fail decisions TCanvas fCanvas; // Canvas for displaying histograms bool fDrawHistos; // Draw histograms? bool fWriteHistos; // Write histograms out? Note that all that members start with a lowercase f, a convention of ours that help recognise data members when reading code.
   
   Classes starting with the letter T (TH1F, TNtuple and TCanvas) are ROOT classes so fQtotHisto, fNhitHisto, fEvtQHisto, fEvtTHisto and fCanvas are all ROOT objects. The ROOT Crib has some introductory information about some ROOT classes that you will need to know about.
   
   Its worth looking at the UserAnalysis constructor UserAnalysis::UserAnalysis() in UserAnalysis.cxx which starts with an initialiser list:- fQtotHisto("fQtotHisto", "Total Charge", 50, 0.0, 1000.0), fNhitHisto("fNhitHisto", "N?Hit!", 50, 0.0, 500.0), fEvtQHisto("fEvtQHisto", "Event Charge", 40, 0.0, 20.0), fEvtTHisto("fEvtTHisto", "Event Times", 100, 0.0, 1000.0), // Ntuple fNtuple("fNtuple","Event Summary",gsChTags), // Module data fNevents(0), fSumQ(0.0), fSumQsqr(0.0), fQminCut(0.0), fQtotMaxCut(100.0), // Classes for drawing fCanvas(std::tmpnam(0),"User Analysis",600,300), fDrawHistos(true), fWriteHistos(false) to ensure that the state is properly initialised. If you are unfamiliar with the initialiser list it might be worth looks at the class section of the C++ Crib

2. Next let's understand what makes UserAnalysis a module, or more strictly a JobModule. From UserAnalysis.h:- class UserAnalysis : public JobCModule { .. }; which means that UserAnalysis inherits from, and hence is a type of JobCModule. In UserAnalysis.cxx you will see the line:- JOBMODULE(UserAnalysis,"UserAnalysis","Compute and record total charge\n"); JOBMODULE is a C++ macro (part of the language, nothing to do with ROOT) that's defined in JobControl/JobCModuleRegistry.h that generates a bit of C++. Its not really necessary to understand what it does, just simply that it will register UserAnalysis with the name "UserAnalysis" and the short description "Compute and record total charge\n".
   
   Now, as UserAnalysis is a JobCModule it has to behave like a JobCModule and respond to all the messages that JobCModule can be sent. If you look at JobControl/JobCModule.h you will see that there are a lot! To simplify matters JobCModule provides defaults, which tell the user that the method is not available so UserAnalysis only has to declare those that it will implement which is does in UserAnalysis.h:- public: UserAnalysis(); // Class constructor // JobCModule methods implemented by this module JobCResult Ana(const MomNavigator *mom); // Analysis method void EndJob(); // End of job method const Registry& DefaultConfig() const; // Default configuration void Config(const Registry& r); // Configure given registry void HandleCommand(JobCommand* cmd); // Handle a job command // Deprecated -- left in for backwards compatibility void HandleCommand(JobCommand* cmd); // Handle a job command .
   
   
3. The DefaultConfig (UserAnalysis::DefaultConfig in UserAnalysis.cxx) method sets up a Registry object with default values:- static Registry r; // Set name of config std::string name = this->GetName(); name += ".config.default"; r.SetName(name.c_str()); // Set values of config r.UnLockValues(); r.Set("QminCut", 0.0); r.Set("QtotCut", 400.0); r.Set("DrawHistos", true); r.Set("WriteHistos",true); r.Set("NbinT", 200); r.Set("Tmin", 0.0); r.Set("Tmax", 200.0E-9); r.LockValues(); return r;
4. The method Config (UserAnalysis::Config in UserAnalysis.cxx) uses a Registry (either the default or one that has been modified at execution time) to change state:- char tmpb; double tmpd; int tmpi; if (r.Get("DrawHistos", tmpb)) { fDrawHistos = tmpb; } if (r.Get("WriteHistos",tmpb)) { fWriteHistos = tmpb; } if (r.Get("QminCut", tmpd)) { fQminCut = tmpd; } if (r.Get("QtotCut", tmpd)) { fQtotCut = tmpd; } ...
5. Next look at the HandleCommand method (UserAnalysis::HandleCommand in UserAnalysis.cxx). Its deprecated because the standard way to configure is via a Registry object, but for a while you may come across this method so its worth looking at. It starts:- if (cmd->HaveCmd()) { // If we have a command... string sc = cmd->PopCmd(); // Get the command if (sc == "Draw") { // If the command is "Draw" fDrawHistos = true; // Default is on if (cmd->HaveOpt()) { // If we have options string sopt = cmd->PopOpt(); // Get the option // If the option is "Off" turn drawing off if (sopt=="Off" || sopt=="off" || sopt=="OFF") fDrawHistos = false; } return; Its uses the JobCommand object it has been passed which contains the module configuration command from the .C file and changes state accordingly
   
   
6. .Now look at the EndJob method (UserAnalysis::EndJob in UserAnalysis.cxx) which computes some totals and writes out all the histograms if required. double aveQ = fSumQ/fNevents; double sigmaQ = sqrt(fSumQsqr/fNevents - aveQ*aveQ); MSG("User",Msg::kInfo) << "Total Events: " << fNevents << "\n"; MSG("User",Msg::kInfo) << " Average Q=" << aveQ << " SigmaQ=" << sigmaQ << "\n"; // Create a file to write the histogram to. By default it becomes // the current directory. Write uses current directory by default if (fWriteHistos) { TFile userFile("useranalysis.root","RECREATE"); fQtotHisto.Write(); fNhitHisto.Write(); fNtuple.Write();
   
   
7. The Ana method (UserAnalysis::Ana in UserAnalysis.cxx) is where all the real work gets done. It gets called for each event and passes a MomNavigator which gives it access to the data. It starts:- const CandDigitListHandle *canddigit = this->GetDigitList(mom); Mark, who wrote this code, likes to explicitly show when an object calls itself. The above line could have been written:- const CandDigitListHandle *canddigit = GetDigitList(mom); but using this, which is always a pointer to the current object makes it crystal clear that this is a member function call and not a global function.
   
   Take a look at UserAnalysis::GetDigitList, where you will see MomNavigator being asked for a CandRecord and CandRecord being asked for a CandDigitListHandle (a list of handles to CandDigits):- CandRecord* candrec = dynamic_cast<CandRecord*> (mom->GetFragment("CandRecord", "PrimaryCandidateRecord")); ... CandDigitListHandle* cdlh = dynamic_cast<CandDigitListHandle*> (candrec->FindCandHandle("CandDigitListHandle", "canddigitlist")); Returning to Ana:, having called GetDigitList, the basic structure is:- if (canddigit) { ... // Accumulate statistics ... // Fill the ntuple ... // If graphics are enabled, draw the histograms ... // Make a filter decision on this event. if (ntuple.qtot<fQtotMaxCut) return JobCResult::kPassed; // Passed cut else return JobCResult::kFailed; // Failed cut } else { // Couldn't find a digit list... MSG("User",Msg::kError) << "No CandDigitList\n"; return JobCResult::kError; } // If we get here somehow, don't make any pass/fail decision on the event return JobCResult::kAOK; } Ana has to return a JobCResult which, for purposes of demonstration will depend on the total charge found in the event. If it is less than the data member fQtotMaxCut then it returns a JobCResult representing a success otherwise a failure. If it cannot find a CandDigitList at all then it returns with an error JobCResult.
   
   Now lets look in more detail to the way it gets the total charge for the event:- double qtot; // Total charge ... qtot = this->GetTotalCharge(canddigit,&nhit); again its sends a message to itself, this time calling GetTotalCharge and passing canddigit, the CandDigitListHandle. Its starts by defining some local variables and reseting some histograms:- int nhit = 0; // Number of hit strips double qtot = 0.0; // Total charge of this event in p.e. fEvtTHisto.Reset(); // Reset event histograms fEvtQHisto.Reset(); // Reset event histograms There then follows the main loop over digits:- // The loop is done with a class that iterates over the list CandDigitHandleItr cdhItr(cdlh->GetDaughterIterator()); for (; cdhItr.IsValid(); cdhItr.Next()) { ... } CandDigitHandleItr is a Navigation iterator that returns a CandDigitHandle pointer to the set its created with, in this case its the daughter list of CandDigitList i.e. the list of CandDigitHandles.
   
   Each CandDigitHandle is processed to get its time and charge and histograms updated if necessary. double t; // Time of single digit double q; // Charge of single digit t = (*cdhItr)->GetTime(); // time stored in units of ns q = (*cdhItr)->GetCharge()/100.0; // charge in pe's // If this digit passes cuts then include it in tallies if (q > fQminCut) { fEvtTHisto.Fill(t,q); // Fill histograms (time weighted by charge) fEvtQHisto.Fill(q); // Fill histograms qtot += q; // Accumulate total charge ++nhit; // Accumulate total number of hits } Finally, returning to Ana, the code to accumulate statistics and update the n-tuple is:- // Accumulate statistics fQtotHisto.Fill(qtot); fNhitHisto.Fill(nhit); fSumQ += qtot; fSumQsqr += qtot*qtot; ++fNevents; // Fill the ntuple ntuple.qtot = qtot; ntuple.nhit = nhit; fNtuple.Fill(&ntuple.qtot);

1. To play safe, in case anything has changed since you made the test release, start by:- cd <myname>/MyTest srt_setup -a
2. Confirm that the build processing is working properly by:- touch Demo/UserAnalysis.h gmake Demo.all "touch"ing UserAnalysis.h, which just changes its modification date to the present time, is sufficient to force srt to rebuild anything that depends, directly or indirectly on that file. It should only take a few moments to rebuild the Demo library libdem.so

3. Having confirmed that the build process is working you are ready to start changing code and rebuilding. Take a look at

      Demo/UserAnalysis.h
      Demo/UserAnalysis.cxx 
   

   to see that you basically understand it, and then here are some things you could put some printout statements to see what is going on, for example in UserAnalysis::GetTotalCharge, after:- t = (*cdhItr)->GetTime(); q = (*cdhItr)->GetCharge()/100.0; add:- cout << "t is " << t << " and q is " << q << endl; and then (from your test directory) rebuild and run:- gmake Demo.all loon -q $MINOS_TUTORIAL_MACROS/write_root_file.C $MINOS_TUTORIAL_DATA/F00018143_0000.mdaq.root Note that the times are in seconds (despite any comment you may see to the contrary), and a more natural unit is nanosecs. Rather than multiply by 1.e-9, it is better not to assume any units as default but instead use Conventions/Munits.h and write:- cout << "t is " << t/Munits::nanosecond << " and q is " << q << endl; If it seems odd to divide by the conversion factor, just think "How many nanoseconds are there in t?".

Contact: Nick West <n.west1@physics.ox.ac.uk>

Security, Privacy, Legal