#ifndef __JDAQTIMESLICE__ #define __JDAQTIMESLICE__ #include #include #include #include #include #include "Jeep/JPrint.hh" #include "JDAQ/JDAQRoot.hh" #include "JDAQ/JDAQPreamble.hh" #include "JDAQ/JDAQSuperFrame.hh" #include "JDAQ/JDAQTimesliceHeader.hh" #include "JDAQ/JDAQEvent.hh" #include "JDAQ/JDAQException.hh" #include "JIO/JSerialisable.hh" #include "JIO/JIOLibrary.hh" namespace KM3NETDAQ { namespace { using JIO::JReader; using JIO::JWriter; } /** * Data time slice. */ class JDAQTimeslice : public JDAQPreamble, public JDAQTimesliceHeader, public std::vector { public: /** * Default constructor. */ JDAQTimeslice() : JDAQPreamble(JLANG::JType()), JDAQTimesliceHeader(), std::vector() {} /** * Constructor. * * \param event DAQ event * \param snapshot use shapshot hits */ JDAQTimeslice(const JDAQEvent& event, const bool snapshot = false) : JDAQPreamble(JLANG::JType()), JDAQTimesliceHeader(event.getDAQChronometer()), std::vector() { using namespace std; map > buffer; if (snapshot) { for (vector::const_iterator hit = event.begin (); hit != event.end (); ++hit) buffer[hit->getModuleID()].push_back(*hit); } else { for (vector::const_iterator hit = event.begin(); hit != event.end(); ++hit) buffer[hit->getModuleID()].push_back(*hit); } for (map >::iterator entry = buffer.begin(); entry != buffer.end(); ++entry) { if (!entry->second.empty()) { sort(entry->second.begin(), entry->second.end()); this->push_back(JDAQSuperFrame(JDAQSuperFrameHeader(getDAQChronometer(), entry->first))); this->rbegin()->add(entry->second.size(), entry->second.data()); } } } /** * Virtual destructor. */ virtual ~JDAQTimeslice() { clear(); } /** * Clear data. */ void clear() { for (iterator i = this->begin(); i != this->end(); ++i) i->clear(); } /** * Assignment operator. * * \param timeslice timeslice * \return this timeslice */ JDAQTimeslice& operator=(const JDAQTimeslice& timeslice) { clear(); setDAQChronometer(timeslice.getDAQChronometer()); for (const_iterator i = timeslice.begin(); i != timeslice.end(); ++i) push_back(*i); return *this; } /** * Add another timeslice. * * \param timeslice timeslice * \return this timeslice */ JDAQTimeslice& add(const JDAQTimeslice& timeslice) { using namespace std; map buffer; for (const_iterator i = this->begin(); i != this->end(); ++i) buffer[i->getModuleIdentifier()] = distance(static_cast(*this).begin(),i); for (JDAQTimeslice::const_iterator i = timeslice.begin(); i != timeslice.end(); ++i) { map::const_iterator p = buffer.find(i->getModuleIdentifier()); if (p != buffer.end()) { JDAQSuperFrame& frame = this->at(p->second); frame.add(*i); sort(frame.begin(), frame.end()); } else this->push_back(*i); } return *this; } /** * Read DAQ timeslice from input. * * \param in JReader * \param timeslice JDAQTimeslice * \return JReader */ friend inline JReader& operator>>(JReader& in, JDAQTimeslice& timeslice) { timeslice.clear(); in >> static_cast (timeslice); in >> static_cast (timeslice); in >> static_cast&>(timeslice); return in; } /** * Write DAQ timeslice to output. * * \param out JWriter * \param timeslice JDAQTimeslice * \return JWriter */ friend inline JWriter& operator<<(JWriter& out, const JDAQTimeslice& timeslice) { out << static_cast (timeslice); out << static_cast (timeslice); out << static_cast&>(timeslice); return out; } /** * Print DAQ Timeslice. * * \param out output stream * \param lpr long print * \return output stream */ ostream& print(ostream& out, const bool lpr = false) const { using namespace std; out << this->ClassName() << newline; out << dynamic_cast (*this) << newline; out << dynamic_cast(*this) << newline; for (JDAQTimeslice::const_iterator frame = this->begin(); frame != this->end(); ++frame) { out << ' ' << setw(5) << frame->getModuleID(); out << ' ' << setw(6) << frame->getLength(); out << ' ' << setw(6) << frame->getDataType(); out << ' ' << setw(6) << frame->size(); if (!lpr) { if (!frame->empty()) { out << ' ' << setw(10) << frame-> begin()->getT(); out << " ... "; out << ' ' << setw(10) << frame->rbegin()->getT(); } out << newline; } else { out << newline; int n = 1; for (JDAQFrame::const_iterator hit = frame->begin(); hit != frame->end(); ++hit, ++n) out << setw(2) << (int) hit->getPMT() << ' ' << setw(8) << hit->getT() << (n%10 == 0 ? newline : whitespace); out << newline; } } return out; } /** * Get size of object. * * \return number of bytes */ virtual int getSize() const { int len = 0; len += JDAQPreamble ::sizeOf(); len += JDAQTimesliceHeader::sizeOf(); len += sizeof(int); for (const_iterator frame = begin(); frame != end(); ++frame) len += frame->getSize(); return len; } ClassDef(JDAQTimeslice,1); }; /** * Print DAQ Timeslice. * * \param out output stream * \param timeslice timeslice * \return output stream */ inline ostream& operator<<(ostream& out, const JDAQTimeslice& timeslice) { return timeslice.print(out, getLongprint(out)); } } #endif