Menu
▾
▴
[Nsclspectcl-cvscommits] SF.net SVN: nsclspectcl:[1721] branches/LLNLMadcChainSpecTcl
|
From: <ro...@us...> - 2010-10-15 11:11:10
|
Revision: 1721
http://nsclspectcl.svn.sourceforge.net/nsclspectcl/?rev=1721&view=rev
Author: ron-fox
Date: 2010-10-15 11:11:04 +0000 (Fri, 15 Oct 2010)
Log Message:
-----------
First reasonable results from multiplicity software.
Modified Paths:
--------------
branches/LLNLMadcChainSpecTcl/CMADC32Unpacker.cpp
branches/LLNLMadcChainSpecTcl/CMADC32Unpacker.h
branches/LLNLMadcChainSpecTcl/Makefile
branches/LLNLMadcChainSpecTcl/MySpecTclApp.cpp
Added Paths:
-----------
branches/LLNLMadcChainSpecTcl/CMultiplicity.cpp
branches/LLNLMadcChainSpecTcl/CMultiplicity.h
Modified: branches/LLNLMadcChainSpecTcl/CMADC32Unpacker.cpp
===================================================================
--- branches/LLNLMadcChainSpecTcl/CMADC32Unpacker.cpp 2010-09-26 13:48:49 UTC (rev 1720)
+++ branches/LLNLMadcChainSpecTcl/CMADC32Unpacker.cpp 2010-10-15 11:11:04 UTC (rev 1721)
@@ -15,6 +15,7 @@
*/
#include <config.h>
#include "CMADC32Unpacker.h"
+#include <TreeParameter.h>
#include <Event.h>
#include <stdint.h>
@@ -61,10 +62,12 @@
//
/*!
- Constrution is a no-op.
+ Constrution just initializes the tree parameter that lets us find the timestamp.
*/
-CMADC32Unpacker::CMADC32Unpacker()
+CMADC32Unpacker::CMADC32Unpacker() :
+ m_timestamp(*(new CTreeParameter("timestamp")))
{
+
}
/*!
@@ -151,18 +154,29 @@
longsRead++;
offset += 2;
}
- // The datum should be the trailer.. verify this.. If so,
- // then save the count field ans parameter 32.
+ // The timestamp should be in the trailer. If so then
+ // set it. There will only be a single tiemstamp per event as we ensure that we
+ // only analyze one module/event in each superevent by breaking out at the first non-event
+ // Events that don't have a trailer won't get a timestamp.
uint32_t timestamp;
if (((datum & ALL_TYPEMASK) >> ALL_TYPESHFT) == TYPE_TRAILER) {
timestamp = datum & TRAILER_COUNTMASK;
+ m_timestamp = timestamp;
}
else {
// Sometimes I see wonky events...missing trailers... data with the top 1/2 filled with only 1's.
- // I'm going to pop an exeption out so that the caller can move on to the next event/supervevent
+ // I'm going to pop an exeption out so that the caller can move on to the next event/supervevent,
+ // If this is not just another header. If it's another header I'll assume the trailer is just
+ // missing.
- throw string("MADC Unpacker did not see a trailer at the end of an MADC32 event");
+ if ((datum & ALL_TYPEMASK >> ALL_TYPESHFT) == TYPE_HEADER) {
+ return offset - 2; // Put the header back in the event for then ext call...
+ }
+ else {
+
+ throw string("MADC Unpacker did not see a trailer at the end of an MADC32 event");
+ }
}
return offset;
}
Modified: branches/LLNLMadcChainSpecTcl/CMADC32Unpacker.h
===================================================================
--- branches/LLNLMadcChainSpecTcl/CMADC32Unpacker.h 2010-09-26 13:48:49 UTC (rev 1720)
+++ branches/LLNLMadcChainSpecTcl/CMADC32Unpacker.h 2010-10-15 11:11:04 UTC (rev 1721)
@@ -21,7 +21,7 @@
#include "CModuleUnpacker.h"
#endif
-
+class CTreeParameter;
/*!
This module is responsible for unpacking the Mesytec MADC-32 32 channel ADC.
These modules have a virtual slot number. It is therefore possible to
@@ -37,6 +37,8 @@
*/
class CMADC32Unpacker : public CModuleUnpacker
{
+private:
+ CTreeParameter& m_timestamp; // Parameter holding the timestamp.
// canonicals:
public:
CMADC32Unpacker();
Added: branches/LLNLMadcChainSpecTcl/CMultiplicity.cpp
===================================================================
--- branches/LLNLMadcChainSpecTcl/CMultiplicity.cpp (rev 0)
+++ branches/LLNLMadcChainSpecTcl/CMultiplicity.cpp 2010-10-15 11:11:04 UTC (rev 1721)
@@ -0,0 +1,665 @@
+/*
+ This software is Copyright by the Board of Trustees of Michigan
+ State University (c) Copyright 2005.
+
+ You may use this software under the terms of the GNU public license
+ (GPL). The terms of this license are described at:
+
+ http://www.gnu.org/licenses/gpl.txt
+
+ Author:
+ Ron Fox
+ NSCL
+ Michigan State University
+ East Lansing, MI 48824-1321
+*/
+
+#include <config.h>
+#include "CMultiplicity.h"
+#include <SpecTcl.h>
+#include <TCLInterpreter.h>
+#include <fstream>
+#include <string>
+#include <stdlib.h>
+#include <unistd.h>
+#include <time.h>
+#include <sys/stat.h>
+#include <sys/types.h>
+
+#include <fstream>
+#include <stdio.h>
+
+using namespace std;
+
+/*!
+ Construction builds the tree parameters and what we can of the CTCLVariable.
+*/
+CMultiplicity::CMultiplicity() :
+ m_bins("multiplicity.bins", 512, "bins"),
+ m_minTime("multiplicity.mintime", 100, "ns"),
+ m_frequency("multiplicity.frequency", 16.0e3, "KHz"),
+ m_cycleCount("multiplicity.cycleCount", 100, "cycles"),
+ m_ROIS("multiplicity.ROIS", false),
+ m_timestamp("timestamp"),
+ m_tZeroSet(false),
+ m_tLast(0)
+{
+ CTreeParameter::BindParameters(); // Make sure timestamp is bound (maybe I don't need this)?
+}
+
+/*!
+ When we are attached to SpecTcl we need to ensure that
+ 1. The Tree variables get bound so that
+ 2. We can load them from the saved file if it exists.
+ 3. We can establish a write trace to maintain the saved file on changes.
+
+ @param rAnalyzer references the TCL analyzer that is associated with the anaklysis case.
+ @return Bool_t
+ @retval kfTRUE - Success.
+ @retval kfFALSE - Some error occured.
+*/
+Bool_t
+CMultiplicity::OnAttach(CAnalyzer& rAnalyzer)
+{
+ try {
+ // Bind my tree variables so I can do what follows:
+
+ SpecTcl* pApi = SpecTcl::getInstance();
+ CTCLInterpreter* pInterp = pApi->getInterpreter();
+ CTreeVariable::BindVariables(*pInterp);
+
+ m_ROIS.Bind(*pInterp);
+
+ readTreeVariables(*pInterp);
+
+ establishTraces(*pInterp);
+
+ return kfTRUE;
+ }
+ catch(...) {
+ std::cerr << "CMultiplicity::OnAttach exception caught\n";
+ return kfFALSE; // Report errors as attach failure for now.
+ }
+}
+
+/*!
+ - Initialize the state system so that we are awaiting the t=0 timestamp.
+ - Pre-compute the ticks/bin
+ - Pre-compute the ticks per cycle.
+ - Obtain the gate object pointers... yelling about any that don't really exist.
+
+Parametrs are ignored but:
+ @return Bool_t
+ @retval kfTRUE - This unconditionally works.
+
+*/
+Bool_t
+CMultiplicity:: OnBegin(CAnalyzer& rAnalyzer,
+ CBufferDecoder& rDecoder)
+{
+ m_tZeroSet = false;
+
+ m_ticksPerBin = (m_bins * m_frequency)/1.0e6;
+ m_ticksPerCycle = m_ticksPerBin * m_bins;
+
+ SpecTcl* pApi = SpecTcl::getInstance();
+
+ bool done = false;
+ int i = 0;
+ while(!done) {
+ char index[100];
+ sprintf(index, "%d", i);
+ const char* pGateName = m_ROIS.Get(TCL_GLOBAL_ONLY, index);
+ if (pGateName) {
+ CGateContainer* pGate = pApi->FindGate(string(pGateName));
+ if (!pGate) {
+ cerr << "Unable to locate ROI named " << pGateName << " Ignored for multiplicity calculations\n";
+ }
+ else {
+ m_Gates.push_back(pGate);
+ }
+ }
+ else {
+ done = true; // no more ROI names.
+ }
+ i++;
+ }
+ return kfTRUE;
+}
+
+/*!
+ Handle end of run buffers. In this case we need to
+ - Write out the multiplicity data in the form required by Phil,
+ - release the resources allocated by the initialization/accumulatioon
+ of these histograms.
+
+ @param rAnalyzer - ignored.
+ @param rBuffer - reference to the buffer decoder. This is used to generate the filename.
+
+ @return Bool_t
+ @retval kfTRUE unconditionally as there's no reason I can think of we should
+ end the pipeline if we have errors.
+*/
+Bool_t
+CMultiplicity::OnEnd(CAnalyzer& rAnalyzer,
+ CBufferDecoder& rBuffer)
+{
+
+ // Get the output filename.
+
+ unsigned int run = rBuffer.getRun();
+ string title = rBuffer.getTitle();
+ string filename = multiplicityFilename(run);
+
+ // Write the data
+
+ writeMultiplicityData(filename, title);
+
+ // Cleanup for next time around:
+
+ m_Gates.clear(); // Need to repopulate this on the next begin run.
+ for (int i = 0; i < m_Multiplicities.size(); i++) {
+ free(m_Multiplicities[i]);
+ }
+ m_Multiplicities.clear();
+
+}
+
+/*!
+ Process other buffer types... we are going to process a special user define buffer
+ that indicates skips in time.
+*/
+Bool_t
+CMultiplicity::OnOther(UInt_t nType,
+ CAnalyzer& rAnalyzer,
+ CBufferDecoder& rDecoder)
+{
+ //TODO: The whole buffer skip thingy.
+
+
+}
+
+/*!
+ Process an Event; If no events have been seen, the base timestamp must be set.
+ The event structure is such that there will always be a timestamp for an event.
+ The gates in m_Gates are evaluated and the trues are counted to give the number of neutrons.
+ The timestamp is converted into a cycle number and bin, and the appropriate bin in the histogram
+ is incremented appropriately (allocated if needed).
+
+ @param pEvent - Raw event -- ignored.
+ @param rEvent - Unpacked event. This is what we use under the assumption a prior event processor
+ has already unpacked the data.
+
+ All other parameters are ignored.
+
+ @return Bool_t
+ @retval kfTRUE - Success
+ @retval kfFALSE - Failure... this will typically be due to an inability to allocated storage.
+
+*/
+Bool_t
+CMultiplicity::operator()(const Address_t pEvent,
+ CEvent& rEvent,
+ CAnalyzer& rAnalyzer,
+ CBufferDecoder& rDecoder)
+{
+ if (m_timestamp.isValid()) { // In case there's un-stamped or empty event.
+
+ // Figure out the stamp relative to the first event and then calculate the
+ // cycle number/bin.
+
+ double stamp = m_timestamp;
+
+ if(!m_tZeroSet) {
+ m_tZero = stamp;
+ m_tZeroSet = true;
+ }
+ stamp = stamp - m_tZero; // TODO: Deal with timestamp wraps.
+ m_tLast = stamp; // Last timestamp received.
+
+ int cycle = static_cast<int>(stamp/m_ticksPerCycle); // Cycle number.
+ int bin = static_cast<int>(stamp/m_ticksPerBin); // Absolute bin number.
+ bin -= cycle*m_bins; // Bin relative to cycle.
+
+ // The neutron multiplicity is just the number of gates that are true for this
+ // event:
+
+ int multiplicity = 0;
+ for (int i =0; i < m_Gates.size(); i++) {
+ CGateContainer* pGate = m_Gates[i];
+ if ( (*pGate)->inGate(rEvent)) {
+ multiplicity++;
+ }
+
+ }
+ if (cycle >= m_Multiplicities.size()) {
+ expandMultiplicities(cycle);
+ }
+ uint32_t* pCycleData = m_Multiplicities[cycle];
+ pCycleData[bin] += multiplicity;
+
+
+ }
+ return kfTRUE;
+
+}
+//////////////////////////////////////////////////////////////////////////////
+//
+// Utilities:
+
+/**
+ * Save the values of the tree variables in a file.
+ * The variables are written as a Tcl script so this really just involves interpreting that script.
+ * Notes:
+ * - We're going to assume the presence of the HOME environement variable.
+ * - The file is hard coded as $HOME/.cycles
+ * - Errors are thrown as exceptions.
+ */
+ void
+ CMultiplicity::saveTreeVariables()
+{
+ double treeValue;
+ std::string gateName;
+ std::string filename = cycleFilename();
+
+ // Open the file and throw on error (e.g. .cycles protected or HOME not writable??).
+
+ std::fstream cyclesFile(filename.c_str(), std::fstream::out);
+ if (cyclesFile.fail()) {
+ throw "Unable to open ~/.cycles file\n";
+ }
+
+ // Start writing:
+
+ // Tree variables:
+
+ cyclesFile << "set multiplicity.bins " << m_bins.getValue() << std::endl;
+ cyclesFile << "set multiplicity.mintime " << m_minTime.getValue() << std::endl;
+ cyclesFile << "set multiplicity.frequency " << m_frequency.getValue() << std::endl;
+ cyclesFile << "set multiplicity.cycleCount " << m_cycleCount.getValue() << std::endl;
+
+ // The gate names also get written to be complete.
+ // The assumption is that multiplicity.ROIS is an array with contiguous integer indices
+ // e.g. 0,1...last...we're going to keep on iterating on the indices until we dont' get
+ // a value. This allows the user to define an arbitrary number of neutron gates in the case that
+ // the detector array grows or shrinks.
+ //
+ int index = 0; // Starting index.
+ char cIndex[100]; // textual version of index.
+ bool done = false;
+
+ // prior to writing the variables the array is unset so that there are no lingering values in case
+ // the array shrinks:
+
+ cyclesFile << "catch {unset multiplicity.ROIS}\n";
+
+ while(!done) {
+ sprintf(cIndex, "%d", index);
+ const char* value = m_ROIS.Get(TCL_GLOBAL_ONLY, cIndex);
+ if(value) {
+ cyclesFile << "set multiplicity.ROIS(" << cIndex << ") " << value << std::endl;
+ }
+ else {
+ done = true;
+ }
+ index++;
+ }
+ // Close the file:
+
+ cyclesFile.close();
+
+}
+/**
+ * Read the tree variables, and gate names from file.
+ * - If there is no ~/.cycles file, this is a no-op.
+ * - The cycles file is a Tcl script it is just interpreted.
+ * - Any errors in interpreting the cycles file are passed to the caller as exceptions.
+ * @param rInter - Reference to the object encapsulated Tcl interpreter.
+ */
+void
+CMultiplicity::readTreeVariables(CTCLInterpreter& rInterp)
+{
+ std::string cyclesFile = cycleFilename();
+ if (!access(cyclesFile.c_str(), R_OK)) {
+
+ rInterp.EvalFile(cyclesFile);
+ }
+
+}
+
+/**
+ * Establishes the write traces on our variables. A single write trace will do because
+ * the only thing it will do is establish object context and invoke saveTreeVariables.
+ */
+void
+CMultiplicity:: establishTraces(CTCLInterpreter& rInterp)
+{
+ Tcl_Interp* pInterp = rInterp.getInterpreter();
+ int flags = TCL_TRACE_WRITES | TCL_GLOBAL_ONLY;
+ Tcl_TraceVar(pInterp, "multiplicity.bins", flags, CMultiplicity::varModified, this);
+ Tcl_TraceVar(pInterp, "multiplicity.mintime", flags, CMultiplicity::varModified, this);
+ Tcl_TraceVar(pInterp, "multiplicity.frequency", flags, CMultiplicity::varModified, this);
+ Tcl_TraceVar(pInterp, "multiplicity.cyclecount", flags, CMultiplicity::varModified, this);
+ Tcl_TraceVar(pInterp, "multiplicity.ROIS", flags, CMultiplicity::varModified, this);
+}
+
+/**
+ * Called when any write trace on the variable is called. If the input flags are not indicating
+ * The interpreter or trace is being destroyed, the variables will be written to file.
+ * @param clientData - Actually a pointer to a CMultiplicity object.
+ * @param interp - Pointer to the Tcl_Interp of the interpreter in which the set happened.
+ * @param name1 - Name of the variable being modified.
+ * @param name2 - If the variable is an array, index of the variable being modified.
+ * @param flags - Trace flags, where TCL_TRACE_DESTROYED and TCL_INTERP_DESTROYED are important to us.
+ */
+char*
+CMultiplicity::varModified(ClientData clientData,
+ Tcl_Interp* pInterp,
+ const char* name1,
+ const char* name2,
+ int flags)
+
+{
+ try {
+ CMultiplicity* pObject = reinterpret_cast<CMultiplicity*>(clientData);
+ pObject->saveTreeVariables();
+ }
+ catch (...) {
+ return const_cast<char*>("Write to ~/.cycles failed\n");
+ }
+ return reinterpret_cast<char*>(0);
+}
+/**
+ * Compute the name of the cycles . file.
+ * - If $HOME is not definedan exception is thrown.
+ */
+std::string
+CMultiplicity::cycleFilename()
+{
+ std::string filename;
+ const char* pHome;
+
+ // Construct and open the file the filename:
+
+ pHome = getenv("HOME");
+ if (!pHome) {
+ throw "HOME Env. variable does not exist!!!\n";
+ }
+ filename = pHome;
+ filename += "/";
+ filename += ".cycles";
+
+ return filename;
+}
+/**
+ * Called to expand the multiplicities histogram array so that the specified cycle can fit.
+ * @param cycle - The largest cycle the histograms must hold.
+ * Note: for now the histograms are just expanded large enough to hold the specified cycle
+ * rather than expanding ahead. The assumption is that the cycle time is such that,
+ * modulo skips, there will always be quite a few events in each cycle.
+ */
+void
+CMultiplicity::expandMultiplicities(int cycle)
+{
+ while (m_Multiplicities.size() <= cycle) {
+ uint32_t* pCycle = reinterpret_cast<uint32_t*>(calloc(m_bins, sizeof(uint32_t)));
+ m_Multiplicities.push_back(pCycle); // TODO: Handle calloc failure.
+ }
+}
+/**
+ * Called to get the name of the multiplicity file. This will be stored in
+ * $HOME/multiplicitydata/run-%04d.multiplicity
+ * @param run - Run number of the run for which we need the filename.
+ * @return std::string
+ * @retval full, absolute path to the file.
+ */
+string
+CMultiplicity::multiplicityFilename(int run)
+{
+ string result;
+
+ // Figure out the directory path:
+
+ const char* home = getenv("HOME");
+ if (!home) {
+ home = "/tmp"; // This has to exist even if HOME does not.
+ }
+ result = home;
+ result += "/multiplicitydata/";
+
+ // Ensurethe directory exists:
+
+ mkdir(result.c_str(), 0777); // let umask determine the permissions. TODO: handle mkdir failures.
+
+
+ // Figure out the filename, append it to the path and return the whole mess to the
+ // caller.
+
+ char filename[100];
+ sprintf(filename, "run-%04d.multiplicity", run);
+ result += filename;
+
+
+ return result;
+
+}
+/**
+ * Write the multiplicity data to file in the format specified by LLNL
+ * This data consists of three sections, each written separately:
+ * - Header that describes the run conditions for the most part.
+ * - Summary that descdribes the average count rates for various chunks of information.
+ * - Counter data for each time width in a cycle.
+ * - Comment section (we'll figure out how to put the title in here.
+ * @param filename - Name of the file to write.
+ * @param comment - Comment string written at the bottom of the file.
+ */
+void
+CMultiplicity::writeMultiplicityData(string filename, string comment)
+{
+
+ std::fstream f(filename.c_str(), ios_base::out); // TODO: handle open failure.
+
+ writeHeader(f, filename);
+ writeSummary(f);
+ for(int i = 0; i < m_bins; i++) {
+ writeGate(f, i);
+ }
+ writeComment(f, comment);
+}
+
+/**
+ * Write the header of the multiplicity file. This consists of the data and time the
+ * data are being written, the file in whih the data are being written.
+ * The time/date used will be 'now'.
+ * @param f - fstream& that is the open file being written to.
+ * @param filename - Full path to the filename being written.
+ */
+void
+CMultiplicity::writeHeader(std::fstream& f, string filename)
+{
+ // Figure out the date/time string:
+
+ char timeStamp[100];
+ struct tm brokenDownTime;
+ time_t epochTime;
+
+ time(&epochTime);
+ localtime_r(&epochTime, &brokenDownTime);
+
+ strftime(timeStamp, sizeof(timeStamp),
+ "%D, %R:%S, data stored",
+ &brokenDownTime);
+
+ f << timeStamp <<endl;
+ f << filename << endl <<endl;
+
+
+}
+/**
+ * Write the summary part of the file header.
+ * - Mode line which is constant.
+ * - Start delay which is the number of microseconds represented by the'
+ * time of the first event.
+ * - Width of bins.
+ * - Number of bins in a cycle
+ * - Total number of cycles memorized.
+ * - Total run time in seconds
+ * - Average count rate for a full cycle.
+ * - Average count rate for the first 20 bins in a cycle.
+ * @param f - fstream& that is open on the output file.
+ */
+void
+CMultiplicity::writeSummary(std::fstream& out)
+{
+ double startDelayms = m_tZero / m_frequency;
+ double runLengthms = (m_tLast - m_tZero) / m_frequency;
+
+ // Now we can write the stuff that does not require looking through the multiplicity data:
+
+ out << "**** Counter Setup ****\n";
+ out << "mode = Neutron Generator-2A\t(Feynman Variance Analysis)\n";
+ out << endl;
+
+ out << startDelayms/1000.0 << ",\t\tstart delay (microseconds)\n";
+ out << ((double)m_minTime/1000.0) << ",\t\tbin width (microseconds)\n";
+ out << (int)((double)m_bins) << ",\t\tbins\n";
+ out << m_Multiplicities.size() << ",\t\ttotal number of cyclees\n";
+ out << endl;
+
+
+ double runLengthSec = runLengthms/1000.0;
+ double binWidthSec = m_minTime / 1.0e9;
+
+
+ // Incorrect text but correct computation:
+
+ out << "Total Run Time = (Cycle Duration) X (Bin Width) X (Total Cycle Count) = Time in seconds\n";
+ out << "Total Run Time = ("
+ << m_bins << ") X "
+ << binWidthSec << " X "
+ << m_Multiplicities.size() << " = "
+ << m_bins * binWidthSec * m_Multiplicities.size() << " seconds\n";
+ out << "\n\n";
+
+ out << "Average count rate for the largest gate ("
+ << (m_bins*m_minTime)*1000 << " microsseconds)\n\n";
+
+ // Need to sum total counts and counts in first 20 bins
+ // There's an implicit assumption that there are at least 20 bins in a cycle.
+
+ int totalCounts = 0;
+ int total20 = 0;
+ for (int i =0; i < m_Multiplicities.size(); i++) {
+ uint32_t* pCycle = m_Multiplicities[i];
+
+ // first 20bins
+ int j = 0;
+ for(; j < 20; j++) {
+ totalCounts += *pCycle;
+ total20 += *pCycle;
+ pCycle++;
+ }
+ for (; j < (int)m_bins; j++) {
+ totalCounts += *pCycle;
+ pCycle++;
+ }
+ }
+ double averageRate = (double)totalCounts/runLengthSec;
+ double averageRate20 = (double)total20/runLengthSec;
+
+ out << "Average Count Rate = Total Count (gate " << (int)((double)m_bins)
+ << " / Total Run Time (sec)\n";
+ out << "Average Count Rate = " << totalCounts << " / "
+ << (double)runLengthSec << " = " << averageRate << endl;
+
+ out << "Average Count Rate for the twentieth gate ( "
+ << (double)m_minTime * 20.0 * 1000.0 << " microseconds)\n\n";
+
+ out << "Average Count Rate = Total Count (gate 20) / Total Run Time (sec)\n";
+ out << "Average Count Rate = " << total20 << " / "
+ << (double)runLengthSec << " = " << averageRate20 << endl <<endl;
+
+}
+/**
+ * Write the data for a gate to the output file. Each gate consists of
+ * a short header indicating how many multiplets we are writing, the gate number
+ * the gate-width (microseconds)
+ * this is followed by a table of multiplet number and occurence counts..
+ * as well as the number of counts seen in that gate over the entire run.
+ * @param out - Reference to the fstream on which the data are to be written
+ * @param gate - Number of bins in the gate.
+ */
+void
+CMultiplicity::writeGate(std::fstream& out, int gate)
+{
+ float gateWidthUsec = m_minTime * gate * 1000.0; // minTime is nsec but gateWidth is in microseconds.
+
+ vector<int> cycleHistogram;
+ int totalCounts = 0;
+
+
+ // Compute the multiplicity histogram for this gate:
+
+ for (int i =0; i < m_Multiplicities.size(); i++) {
+ uint32_t* pCycle = m_Multiplicities[i];
+
+ int gateSum = 0;
+ for (int j = 0; j < gate; j++) {
+ totalCounts += *pCycle;
+ gateSum += *pCycle++;
+
+ }
+ if (cycleHistogram.size() <= gateSum) expandHistogram(cycleHistogram, gateSum);
+ cycleHistogram[gateSum]++;
+ }
+ int multipletCount = cycleHistogram.size();
+ if (multipletCount < 40) multipletCount = 40;
+
+ // Write the header:
+
+ out << "***** Counter Data *****\n\n";
+ out << multipletCount << ",\t\tTotal number of multplet entries\n";
+ out << gate << ",\t\tgate\n";
+ out << gateWidthUsec << ",\t\tgate-width (microseconds)\n\n";
+
+
+ // write the data
+
+ out << "Multiplet No.,\tData\n\n";
+ for (int i =0; i < multipletCount; i++) {
+ int value = 0;
+ if (i < cycleHistogram.size()) value = cycleHistogram[i];
+
+ out << i << ",\t\t" << value << endl;
+ }
+ out << endl;
+
+ // Write the sum:
+
+ out << "Sum counts = \t" << totalCounts << endl << endl << endl;
+}
+
+/**
+ * Expands a vector so that it can hold a specific max index.
+ * @param histogram - vector<int>& to expand.
+ * @param size - largest index the vector should hold when we are done.
+ * the vector is zero filled.
+ */
+void
+CMultiplicity::expandHistogram(vector<int>& histogram, int size)
+{
+ while(size >= histogram.size()) histogram.push_back(0);
+}
+/**
+ * Write the comment part of the multiplicity output file.
+ * @param f - Reference to fstream open on output file.
+ * @param comment - Comment to write to the comment segment of the output file.
+ */
+void
+CMultiplicity::writeComment(std::fstream& f, string comment)
+{
+ f << endl << endl << endl;
+ f << "*****User comments*****\n";
+ f << "\t" << comment;
+
+}
Added: branches/LLNLMadcChainSpecTcl/CMultiplicity.h
===================================================================
--- branches/LLNLMadcChainSpecTcl/CMultiplicity.h (rev 0)
+++ branches/LLNLMadcChainSpecTcl/CMultiplicity.h 2010-10-15 11:11:04 UTC (rev 1721)
@@ -0,0 +1,145 @@
+/*
+ This software is Copyright by the Board of Trustees of Michigan
+ State University (c) Copyright 2005.
+
+ You may use this software under the terms of the GNU public license
+ (GPL). The terms of this license are described at:
+
+ http://www.gnu.org/licenses/gpl.txt
+
+ Author:
+ Ron Fox
+ NSCL
+ Michigan State University
+ East Lansing, MI 48824-1321
+*/
+
+
+#ifndef __CMULTIPLICITY_H
+#define __CMULTIPLICITY_H
+
+#ifndef __EVENTPROCESSOR_H
+#include <EventProcessor.h>
+#endif
+
+
+#ifndef __TREEPARAMETER_H
+#include <TreeParameter.h>
+#endif
+
+#ifndef __TCLVARIABLE_H
+#include <TCLVariable.h>
+#endif
+
+
+#ifndef __TCL_H
+#include <tcl.h>
+#ifndef __TCL_H
+#define __TCL_H
+#endif
+#endif
+
+
+#ifndef __STL_VECTOR
+#include <vector>
+#ifndef __STL_VECTOR
+#define __STL_VECTOR
+#endif
+#endif
+
+#ifndef __STL_STRING
+#incldue <string>
+#ifndef __STL_STRING
+#define __STL_STRING
+#endif
+#endif
+
+#ifndef __CRTL_STDINT_H
+#include <stdint.h>
+#ifndef __CRTL_STDINT_H
+#define __CRTL_STDINT_H
+#endif
+#endif
+
+
+class CTCLInterpreter;
+class CGateContainer;
+class fstream;
+/*!
+ Event processor to handle the multiplicity calculations.
+ Assumptions:
+ - Data from only one module is present in the event, which implies:
+ - There's only a unique timestamp.
+ Several tree parameters steer this:
+ - 'bins' - (integer, no units) Number of bins in a cycle
+ - 'minTime' - (floating nanoseconds) Minimum time per bin (defines the length of the cycle.
+ - 'frequency' - (floating Khz) Frequency of the timestamp timebase pulser.
+ - 'cycleCount' - (integer) Number of cycles that must be accepted prior to allowing data to skip.
+ - 'ROIs' - Array of 24 strings that are the names of the gates that define what a neutron is.
+*/
+class CMultiplicity : public CEventProcessor
+{
+
+ // Local data:
+private:
+ CTreeVariable m_bins;
+ CTreeVariable m_minTime;
+ CTreeVariable m_frequency;
+ CTreeVariable m_cycleCount;
+ CTCLVariable m_ROIS;
+ CTreeParameter m_timestamp; // The event parameter holding the timestamp.
+ std::vector<CGateContainer*> m_Gates;
+ bool m_tZeroSet; // True if first event seen to set time zero
+ int m_tZero; // Timestamp value at t==0.
+ int m_tLast; // Last timestamp.
+ double m_ticksPerBin; // Floating to avoid potential aliasing in case of non-conformalism
+ double m_ticksPerCycle;
+ std::vector<uint32_t*> m_Multiplicities; // Multiplicity histograms.
+
+
+
+
+ // Canonicals
+public:
+ CMultiplicity();
+
+ // Interface definitions:
+public:
+ virtual Bool_t OnAttach(CAnalyzer& rAnalyzer); // Called on registration.
+ virtual Bool_t OnBegin(CAnalyzer& rAnalyzer,
+ CBufferDecoder& rDecoder); // Begin Run.
+ virtual Bool_t OnEnd(CAnalyzer& rAnalyzer,
+ CBufferDecoder& rBuffer); // End Run.
+ virtual Bool_t OnOther(UInt_t nType,
+ CAnalyzer& rAnalyzer,
+ CBufferDecoder& rDecoder); // Unrecognized buftype.
+
+ virtual Bool_t operator()(const Address_t pEvent,
+ CEvent& rEvent,
+ CAnalyzer& rAnalyzer,
+ CBufferDecoder& rDecoder);
+ // Utitities:
+private:
+ void saveTreeVariables();
+ void readTreeVariables(CTCLInterpreter& rInterp);
+ void establishTraces(CTCLInterpreter& rInterp);
+ void expandMultiplicities(int cycle);
+ void writeHeader(std::fstream& out, std::string filename);
+ void writeSummary(std::fstream& out);
+ void writeMultiplicityData(std::string filename, std::string comment);
+ void writeGate(std::fstream& f, int gate);
+ void writeComment(std::fstream& f, std::string comment);
+
+ static char* varModified(ClientData clientData,
+ Tcl_Interp* pInterp,
+ const char* name1,
+ const char* name2,
+ int flags);
+ static std::string cycleFilename();
+ static std::string multiplicityFilename(int run);
+ static void expandHistogram(std::vector<int>& histogram, int size);
+
+};
+
+
+#endif
Modified: branches/LLNLMadcChainSpecTcl/Makefile
===================================================================
--- branches/LLNLMadcChainSpecTcl/Makefile 2010-09-26 13:48:49 UTC (rev 1720)
+++ branches/LLNLMadcChainSpecTcl/Makefile 2010-10-15 11:11:04 UTC (rev 1721)
@@ -29,7 +29,7 @@
CRateProcessor.o CRateList.o CRateCommand.o CFitButton.o \
CStackMapCommand.o CModuleUnpacker.o C785Unpacker.o CNADC2530Unpacker.o \
CStackUnpacker.o CMADC32Unpacker.o CRateEventProcessor.o \
- CV1x90Unpacker.o CV977Unpacker.o CMASE.o
+ CV1x90Unpacker.o CV977Unpacker.o CMASE.o CMultiplicity.o
#
Modified: branches/LLNLMadcChainSpecTcl/MySpecTclApp.cpp
===================================================================
--- branches/LLNLMadcChainSpecTcl/MySpecTclApp.cpp 2010-09-26 13:48:49 UTC (rev 1720)
+++ branches/LLNLMadcChainSpecTcl/MySpecTclApp.cpp 2010-10-15 11:11:04 UTC (rev 1721)
@@ -33,6 +33,7 @@
#include "CParamMapCommand.h"
#include "CStackMapCommand.h"
#include "CStackUnpacker.h"
+#include "CMultiplicity.h"
#include "CRateCommand.h"
#include "CRateEventProcessor.h"
#include "CFitButton.h"
@@ -94,7 +95,10 @@
CMySpecTclApp::CreateAnalysisPipeline(CAnalyzer& rAnalyzer)
{
RegisterEventProcessor(*(new CStackUnpacker), "adc-data");
+ RegisterEventProcessor(*(new CMultiplicity), "Multiplicity");
RegisterEventProcessor(*(new CScalerProcessor), "Scalers");
+
+
// RegisterEventProcessor(*(new CRateEventProcessor), "rate-stripchart");
}
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|