Index: applgrid/trunk/src/appl_file.cxx =================================================================== --- applgrid/trunk/src/appl_file.cxx (revision 0) +++ applgrid/trunk/src/appl_file.cxx (revision 1895) @@ -0,0 +1,212 @@ +/** + ** @file appl_file.cxx + ** + ** @author sutt + ** @date Wed 30 Dec 2020 15:25:42 GMT + ** + ** $Id: appl_file.cxx, v0.0 Wed 30 Dec 2020 15:25:42 GMT sutt $ + ** + ** Copyright (C) 2020 sutt (sutt@cern.ch) + ** + **/ + + +#include "appl_grid/appl_file.h" +#include + + +appl::file::file( const std::string& name, const std::string& option ) : + mfilename(name), mfile(0), mopen(false), msize(0) { + + mopt = option; + + /// make sure it is always binary mode + if ( mopt.find("b")==std::string::npos ) mopt += "b"; + + // std::cout << "appl::file::opening file: " << filename() << "\toptions: " << mopt << std::endl; + + /// file sizes - actual filesize and uncopmpressed size + int filesize = 0; + int usize = 0; + + if ( mopt.find("r")!=std::string::npos ) { + + /// first get file size ... + + struct stat stat_buf; + int rc = stat(filename().c_str(), &stat_buf); + if ( rc == 0 ) filesize = stat_buf.st_size; + else { + std::cerr << "appl::file: can not determine file size for file: " << filename() << std::endl; + return; + } + + /// next check if it is in zipped format ... + + FILE* tmp_file = fopen( filename().c_str(), "rb" ); + int zip_signature = 0; + fread( &zip_signature, sizeof(int), 1, tmp_file ); + if ( (zip_signature&0xffffff)==0x88b1f ) { + /// if it is read the file size, and the uncompressed filesize from the file ... + /// uncompressed size + int offset = filesize - 4; + fseek( tmp_file, offset, SEEK_SET ); + fread( &usize, sizeof(int), 1, tmp_file ); + } + else usize = filesize; + + // std::cout << "\tfile size: " << filesize << std::endl; + // std::cout << "\tuncompressed size: " << usize << std::endl; + + fclose( tmp_file ); + + } + + /// now open the file ... + + mfile = gzopen( mfilename.c_str(), mopt.c_str() ); + + mopen = true; + + /// write header + if ( mopt.find("w")!=std::string::npos ) { + + double v = SB::WRITEGUARD; + gzwrite( mfile, (void*)&v, sizeof(double) ); + mindex.add( "header", sizeof(double) ); + + msize += sizeof(double); + + /// the index is added at the end before the file is closed + } + if ( mopt.find("r")!=std::string::npos ) { + + /// read the header ... + + double v = 0; + gzread( mfile, (void*)&v, sizeof(double) ); + + /// if this one of our files ? + if ( v!=SB::WRITEGUARD ) { + std::cerr << "appl::file: incorrect file format file: " << filename() << std::endl; + Close(); + mopen = false; + return; + } + + /// yes it is so ... + + /// first get the (uncompressed) file size from the last 4 bytes + /// of the compressed file, then use that get the offsets to the + /// index stored at the end of the file - *why* does this have to + /// be soooo difficult + + // std::cout << "opened file: " << filename() << "\tsize: " << filesize << std::endl; + + /// read the uncompressed file size from the comprfessed file + /// so that we can navigate through the commpressed file, + /// to get the index in case we want unsequential access + + // std::cout << "uncompressed size: " << usize << std::endl; + + /// get the total file size and offset of the index offset + + double vtrailer[3]; + + gzseek( mfile, usize-sizeof(double)*3, SEEK_SET ); + gzread( mfile, (void*)vtrailer, sizeof(double)*3 ); + + int index_size = (vtrailer[1]-vtrailer[0]-sizeof(double)*3); /// in bytes + + std::vector vindex(index_size/sizeof(double)); + + int indexptr = vtrailer[0]; + + /// std::cout << "vtrailer[0]: " << vtrailer[0] << std::endl; + /// std::cout << "vtrailer[1]: " << vtrailer[1] << std::endl; + /// std::cout << "vtrailer[2]: " << vtrailer[2] << std::endl; + + // gzrewind( mfile ); + gzseek( mfile, indexptr, SEEK_SET ); + gzread( mfile, (void*)&vindex[0], index_size ); + + mindex.deserialise( vindex ); + + // std::cout << mindex << std::endl; + + /// rewind to open the file for proper reading ... + gzrewind( mfile ); + /// read the header again, so that we are aligned with the + /// first object in the file + gzread( mfile, (void*)&v, sizeof(double) ); + + // std::cout << "write guard " << v << std::endl; + + } +} + + +appl::file::~file() { Close(); } + + +void appl::file::Close() { + + if ( isopen() ) { + + /// if writing the file, update the global size + + if ( mopt.find("w")!=std::string::npos ) { + double content_trailer = SB::WRITEGUARD; + int bytes = gzwrite( mfile, (void*)&content_trailer, sizeof(double) ); + mindex.add( "trailer", bytes ); + + msize += bytes; + + double index_offset = msize; + + // std::cout << "writing: index offset: " << index_offset << std::endl; + + /// this is a thorny one - should we add the index to the index ? + /// until the index has been written, it should not be in the index + /// but once it has been written, then it *should* be in the index + /// so we would need to add the index to itself, *before* writing + /// to the file + /// basically, I think it is more sensible *not* to include the + /// index itself into the index, but of course, that could could + /// be changed if necessary + /// what we *should* do, is to add some pointer to the index at + /// the start of the file, so that we can navigate to it directly + /// on openning the file + Write( mindex ); + + /// this doesn;t work - need to work out why + // gzrewind( mfile ); + // gzseek( mfile, sizeof(double), SEEK_SET ); + // double tsize = msize; + + double trailer[3] = { 0, 0, 0 }; + trailer[0] = index_offset; + trailer[1] = msize+3*sizeof(double); + trailer[2] = SB::WRITEGUARD; + + bytes = gzwrite( mfile, (void*)&trailer, sizeof(double)*3 ); + // std::cout << "writing trailer: " << std::endl; + // std::cout << "\tindex offset: " << trailer[0] << std::endl; + // std::cout << "\ttotal size: " << trailer[1] << std::endl; + // std::cout << "\twrite guard: " << trailer[2] << std::endl; + + mindex.add("file trailer", bytes ); + } + + gzclose(mfile); + mopen = false; + + mindex.clear(); + + } + +} + + + + Index: applgrid/trunk/src/equals.h =================================================================== --- applgrid/trunk/src/equals.h (revision 0) +++ applgrid/trunk/src/equals.h (revision 1895) @@ -0,0 +1,37 @@ +/* emacs: this is -*- c++ -*- */ +/** + ** @file equals.h + ** + ** @author sutt + ** @date Wed 30 Dec 2020 11:26:54 GMT + ** + ** $Id: equals.h, v0.0 Wed 30 Dec 2020 11:26:54 GMT sutt $ + ** + ** Copyright (C) 2020 sutt (sutt@cern.ch) + ** + **/ + + +#ifndef EQUALS_H +#define EQUALS_H + + +inline bool equals( double a, double b, double limit=1 ) { + if ( a==b ) return true; + if ( std::fabs(a-b)<(std::fabs(std::min(a,b))*std::numeric_limits::epsilon()*limit ) ) return true; + if ( ( a==0 || b==0 ) && ( std::fabs(a-b)::epsilon()*limit ) ) return true; + return false; +} + + +#endif // EQUALS_H + + + + + + + + + + Index: applgrid/trunk/appl_grid/appl_file.h =================================================================== --- applgrid/trunk/appl_grid/appl_file.h (revision 0) +++ applgrid/trunk/appl_grid/appl_file.h (revision 1895) @@ -0,0 +1,181 @@ +/* emacs: this is -*- c++ -*- */ +/** + ** @file appl_file.h + ** + ** @author sutt + ** @date Wed 30 Dec 2020 15:25:42 GMT + ** + ** $Id: appl_file.h, v0.0 Wed 30 Dec 2020 15:25:42 GMT sutt $ + ** + ** Copyright (C) 2020 sutt (sutt@cern.ch) + ** + **/ + + +#ifndef APPL_FILE_H +#define APPL_FILE_H + +#include + +#include "appl_grid/serialise_base.h" +#include "appl_grid/stream_vector.h" +#include "appl_grid/file_index.h" + +// #include + +namespace appl { + +class file { + +public: + + file( const std::string& name, const std::string& option="rb" ); + + ~file(); + + void Close(); + + void close() { Close(); } + + /// accessors + + std::string filename() const { return mfilename; } + + bool isopen() const { return mopen; } + + /// provide the index + + const file_index& index() const { return mindex; } + + /// object writing ... + +#if 0 + template + void Write( const std::vector& v ) { + std::cout << "write: (0) " << t.name() << std::endl; + if ( !isopen() ) return; + double vsize = v.size(); + double Tsize = sizeof(T); + msize += gzwrite( mfile, (void*)&vsize, sizeof(double) ); + msize += gzwrite( mfile, (void*)&Tsize, sizeof(double) ); + msize += gzwrite( mfile, (void*)(&v[0]), sizeof(T)*v.size() ); + } +#endif + + template + void Write( const T& t ) { + if ( !isopen() ) return; + std::vector v = t.serialise(); + size_t bytes = gzwrite( mfile, (void*)(&v[0]), sizeof(double)*v.size() ); + if ( bytes!=sizeof(double)*v.size() ) std::cerr << "issue writing object " << t.name() << std::endl; + msize += bytes; + mindex.add( t.name(), bytes ); + } + + + template + void Write( const std::string& name, const std::vector& t ) { + stream_vector sv( name, t ); + Write( sv ); + } + + /// object reading ... + + template + void Read( std::vector& v ) { + if ( !isopen() ) return; + double vsize = 0; + double Tsize = 0; + gzread( mfile, &vsize, sizeof(double) ); + gzread( mfile, &Tsize, sizeof(double) ); + v.resize(vsize); + gzread( mfile, &v[0], Tsize*vsize ); + } + + + + template + void Read( T& t ) { + if ( !isopen() ) return; + std::vector v(2); + gzread( mfile, &v[0], sizeof(double)*2 ); + v.resize( v[1] ); + int bytes = gzread( mfile, &v[2], sizeof(double)*(v[1]-2) ); + if ( bytes!=(sizeof(double)*(v[1]-2)) ) std::cerr << "issue writing object " << t.name() << std::endl; + t = T(v); + } + + + template + T Read() { + T t; + Read(t); + return t; + } + + template + T Read( const std::string& name ) { + return Get(name); + } + + /// non-sequential access using the index ... + + template + T Get( const std::string& name ) { + if ( isopen() ) { + file_index::entry e = mindex.find( name ); + if ( e.size>0 ) { + // gzrewind( mfile ); + gzseek( mfile, e.offset, SEEK_SET ); + return Read(); + } + } + std::cerr << "WARNING: could not retrieve object: " << name << std::endl; + return T(); + } + + +private: + + /// filename + std::string mfilename; + + /// options, read write, etc + std::string mopt; + + /// the actual file + gzFile mfile; + + /// is is open ? + bool mopen; + + //the + size_t msize; + + /// file index + file_index mindex; + +}; + +} + + +inline std::ostream& operator<<( std::ostream& s, const appl::file& t ) { + s << "file: " << t.filename() << "\t"; + if ( !t.isopen() ) s << " file is not open"; + else s << t.index(); + return s; +} + + +#endif // APPL_FILE_H + + + + + + + + + + Index: applgrid/trunk/appl_grid/histogram.h =================================================================== --- applgrid/trunk/appl_grid/histogram.h (revision 0) +++ applgrid/trunk/appl_grid/histogram.h (revision 1895) @@ -0,0 +1,309 @@ +/** emacs: this is -*- c++ -*- **/ +/** + ** @file histogram.h + ** + ** @author sutt + ** @date Mon 28 Dec 2020 19:33:12 GMT + ** + ** $Id: histogram.h, v0.0 Mon 28 Dec 2020 19:33:12 GMT sutt $ + ** + ** Copyright (C) 2020 sutt (sutt@cern.ch) + ** + **/ + + +#ifndef HISTOGRAM_H +#define HISTOGRAM_H + +#include +#include +#include +#include +#include + + +#include "serialise_base.h" +#include "binary_search.h" + + + +class histogram { + +public: + + // histogram error exception + class exception : public std::exception { + public: + exception(const std::string& s) { std::cerr << what() << " " << s << std::endl; }; + virtual const char* what() const throw() { return ""; } + }; + +public: + + histogram( const std::string& name="" ) : mname(name) { } + + histogram( const std::string& name, size_t nbins, const double* limits ) : mname(name) { + if ( nbins>0 ) create( nbins, limits ); + else throw exception("histogram: not enough bins creating histogram: "+name); + } + + histogram( const std::string& name, const std::vector& limits ) : mname(name) { + if ( limits.size()>1 ) create( limits.size()-1, &limits[0] ); + else throw exception("histogram: not enough bin limits creating histogram: "+name); + } + + histogram( const std::string& name, size_t nbins, double lo, double hi ) : mname(name) { + if ( nbins==0 ) throw exception("histogram: not enough bins creating histogram: "+name); + std::vector limits(nbins+1,0); + for ( size_t i=nbins ; i-- ; ) limits[i] = (lo*(nbins-i) + hi*i)/nbins; + /// make sure the limits at least are exact + limits[0] = lo; + limits[nbins] = hi; + create( nbins, &limits[0] ); + } + + /// constructor from stream + + histogram( const std::vector stream ) { + deserialise( stream ); + } + + virtual ~histogram() { } + + histogram( const histogram& h ) : + mname(h.mname), + mxlimits(h.mxlimits), + mx(h.mx), + my(h.my), + mye(h.mye) + { } + + + histogram& operator=( histogram&& h ) { + std::swap( mname, h.mname ); + + std::swap( mxlimits, h.mxlimits ); + std::swap( mx, h.mx ); + + std::swap( my, h.my ); + std::swap( mye, h.mye ); + + return *this; + } + + + histogram& operator=( const histogram& h ) { + mname = h.mname; + mxlimits = h.mxlimits; + mx = h.mx; + my = h.my; + mye = h.mye; + return *this; + } + + /// accessors + + std::string name() const { return mname; } + + size_t size() const { return mx.size(); } + + double operator[]( int i ) const { return my.at(i); } + double& operator[]( int i ) { return my[i]; } + + double& y( int i ) { return my.at(i); } + double& ye( int i ) { return mye.at(i); } + + double y( int i ) const { return my.at(i); } + double ye( int i ) const { return mye.at(i); } + + double lo( int i ) const { return mxlimits.at(i); } + double hi( int i ) const { return mxlimits.at(i+1); } + + const std::vector& xlimits() const { return mxlimits; } + + /// logical operators - + /// NB: these only test on the BIN LIMITS, NOT THE BIN CONTENTS !! + + bool operator==( const histogram& h ) { + if ( size()!=h.size() ) return false; + if ( size()==0 ) return true; + double delta = mxlimits[1]-mxlimits[0]; + for ( size_t i=mxlimits.size() ; i-- ; ) if ( std::fabs(mxlimits[i]-h.mxlimits[i])>delta*1e-10 ) return false; + return true; + } + + bool operator!=( const histogram& h ) { + return !(operator==(h)); + } + + /// mathematical operations ... + + histogram& operator*=( double d ) { + for ( size_t i=size() ; i-- ; ) { + my[i] *= d; + mye[i] *= d; + } + return *this; + } + + histogram operator*( double d ) const { + histogram h(*this); + return (h*=d); + } + + histogram& operator/=( double d ) { + return operator*=(1/d); + } + + histogram operator/( double d ) const { + histogram h(*this); + return h/=d; + } + + + histogram& operator+=( const histogram& h ) { + if ( ! operator==(h) ) throw exception( "histogram: bin mismatch for operator +- " ); + for ( size_t i=size() ; i-- ; ) { + my[i] += h.my[i]; + mye[i] = std::sqrt( h.mye[i]*h.mye[i] + mye[i]*mye[i] ); + } + return *this; + } + + histogram operator+( histogram h ) const { + return ( h += (*this) ); + } + + + histogram operator-=( const histogram& h ) { + return operator+=( h*(-1) ); + } + + histogram operator-( const histogram& h ) const { + histogram hi(*this); + return hi-=h; + } + + + /// filling ... + + int index( double x ) const { return find_index( x, mxlimits ); } + + + void fill( double x, double w=1 ) { + int ibin = index( x ); + if ( ibin != -1 ) { + my[ibin] += w; + mye[ibin] = std::sqrt( mye[ibin]*mye[ibin] + w ); + } + } + + + /// serialisation and deserialisation ... + + std::vector serialise() const { + + std::vector s; + + s.push_back( SB::WRITEGUARD ); + + s.push_back(0); // overall size + + SB::serialise( s, name() ); + SB::serialise( s, mxlimits ); + /// don't need to save this - can compute from the limits + /// SB::serialise( s, mx ); + SB::serialise( s, my ); + SB::serialise( s, mye ); + + s.push_back( SB::WRITEGUARD ); + + s[1] = s.size(); + + return s; + } + + + void deserialise( const std::vector& v ) { + + std::vector::const_iterator itr = v.begin(); + + if ( (*itr++) != SB::WRITEGUARD ) throw exception("histogram: read error"); + + /// don't need this now ... + // int ntotal = (*itr++); + (*itr++); + + SB::deserialise( itr, mname ); + SB::deserialise( itr, mxlimits ); + /// SB::deserialise( itr, mx ); + SB::deserialise( itr, my ); + SB::deserialise( itr, mye ); + + double var = (*itr++); + + // std::cout << "read: " << mname << "\tsize: " << ntotal << std::endl; + + if ( var != SB::WRITEGUARD ) throw exception("histogram: read error reading histogram: "+name()); + + mx.resize(my.size()); + for ( size_t i=my.size() ; i-- ; ) mx[i] = 0.5*(mxlimits[i]+mxlimits[i+1]); + + } + +protected: + + void create( size_t nbins, const double* limits ) { + + mxlimits.resize(nbins+1); + + for ( size_t i=nbins+1 ; i-- ; ) mxlimits[i] = limits[i]; + + mx.resize(nbins); + + for ( size_t i=nbins ; i-- ; ) mx[i] = 0.5*(limits[i]+limits[i+1]); + + my = std::vector(nbins,0); + mye = std::vector(nbins,0); + + } + +protected: + + std::string mname; + + std::vector mxlimits; + + std::vector mx; + + std::vector my; + std::vector mye; + +}; + + +inline histogram operator*( double d, histogram h ) { + return (h*=d); +} + + + + +inline std::ostream& operator<<( std::ostream& s, const histogram& h ) { + s << "key: " << h.name() << std::endl; + for ( size_t i=0 ; i +#include + +#include "appl_grid/vector_stream.h" +#include "appl_grid/serialise_base.h" + + + +class stream_grid { + +public: + + // grid error exception + class exception : public std::exception { + public: + exception(const std::string& s) { std::cerr << what() << " " << s << std::endl; }; + virtual const char* what() const throw() { return ""; } + }; + +public: + + stream_grid() { } + + stream_grid( const std::string& name, + const std::vector& xv, + const std::vector& yv, + const std::vector& zv ) : + mname(name), + mxaxis( xv ), + myaxis( yv ), + mzaxis( zv ), + mnz(zv.size()), + mnyz(zv.size()*yv.size()) + { + // mnodes = std::vector > >( xv.size(), std::vector >( yv.size(), std::vector(zv.size(), 0 ) ) ); + mnodes = std::vector(xv.size()*yv.size()*zv.size(), 0 ); + } + + stream_grid( const stream_grid& sg ) : + mname(sg.mname), + mxaxis(sg.mxaxis), + myaxis(sg.myaxis), + mzaxis(sg.mzaxis), + mnz(sg.mnz), + mnyz(sg.mnyz), + mnodes(sg.mnodes) + { } + + + stream_grid& operator=( stream_grid&& sg ) { + std::swap(mname, sg.mname); + std::swap(mxaxis, sg.mxaxis); + std::swap(myaxis, sg.myaxis); + std::swap(mzaxis, sg.mzaxis); + std::swap(mnz, sg.mnz); + std::swap(mnyz, sg.mnyz); + std::swap(mnodes, sg.mnodes); + return *this; + } + + + stream_grid( const std::vector& v ) { + deserialise( v ); + } + + virtual ~stream_grid() { } + + std::string name() const { return mname; } + + // void set( int ix, int iy, int iz, double d ) { mnodes[ix+mnx*iy+mnxy*iz] = d; } + void set( int ix, int iy, int iz, double d ) { mnodes[ix*mnyz+iy*mnz+iz] = d; } + + // double& get( int ix, int iy, int iz ) { return mnodes[ix+mnx*iy+mnxy*iz]; } + // double get( int ix, int iy, int iz ) const { return mnodes[ix+mnx*iy+mnxy*iz]; } + + // double& operator()( int ix, int iy, int iz ) { return mnodes[ix+mnx*iy+mnxy*iz]; } + // double operator()( int ix, int iy, int iz ) const { return mnodes[ix+mnx*iy+mnxy*iz]; } + + double& get( int ix, int iy, int iz ) { return mnodes[ix*mnyz+iy*mnz+iz]; } + double get( int ix, int iy, int iz ) const { return mnodes[ix*mnyz+iy*mnz+iz]; } + + double& operator()( int ix, int iy, int iz ) { return mnodes[ix*mnyz+iy*mnz+iz]; } + double operator()( int ix, int iy, int iz ) const { return mnodes[ix*mnyz+iy*mnz+iz]; } + + const std::vector& xaxis() const { return mxaxis; } + const std::vector& yaxis() const { return myaxis; } + const std::vector& zaxis() const { return mzaxis; } + + const double& xaxis(int i) const { return mxaxis[i]; } + const double& yaxis(int i) const { return myaxis[i]; } + const double& zaxis(int i) const { return mzaxis[i]; } + + std::vector serialise() const { + + std::vector s; + + s.push_back( SB::WRITEGUARD ); + + s.push_back(0); // overall size + + SB::serialise( s, name() ); + SB::serialise( s, mxaxis ); + SB::serialise( s, myaxis ); + SB::serialise( s, mzaxis ); + + SB::serialise( s, mnodes ); + + s.push_back( SB::WRITEGUARD ); + + s[1] = s.size(); + + return s; + } + + + void deserialise( const std::vector& v ) { + + std::vector::const_iterator itr = v.begin(); + + if ( (*itr++) != SB::WRITEGUARD ) throw exception("stream_grid: read error"); + + /// size_t ntotal = (*itr++); + + (*itr++); + + SB::deserialise( itr, mname ); + SB::deserialise( itr, mxaxis ); + SB::deserialise( itr, myaxis ); + SB::deserialise( itr, mzaxis ); + SB::deserialise( itr, mnodes ); + + double var = (*itr++); + + mnz = mzaxis.size(); + mnyz = mzaxis.size()*myaxis.size(); + + if ( var != SB::WRITEGUARD ) throw exception("stream_grid: read error reading grid: "+name()); + + /// if ( ntotal!=v.size() ) std::cerr << "WARNING: mismatch stream_grid size: " << ntotal << " " << v.size() << std::endl; + + } + + + +private: + + std::string mname; + + std::vector mxaxis; + std::vector myaxis; + std::vector mzaxis; + + size_t mnz; + size_t mnyz; + + std::vector mnodes; + +}; + + +inline std::ostream& operator<<( std::ostream& s, const stream_grid& sg ) { + + s << "[ x axis: " << sg.xaxis().size() << "\t" << sg.xaxis()[0] << " .. " << sg.xaxis()[sg.xaxis().size()-1] << " ]\n"; + s << "[ y axis: " << sg.yaxis().size() << "\t" << sg.yaxis()[0] << " .. " << sg.yaxis()[sg.yaxis().size()-1] << " ]\n"; + s << "[ z axis: " << sg.zaxis().size() << "\t" << sg.zaxis()[0] << " .. " << sg.zaxis()[sg.zaxis().size()-1] << " ]\n"; + + +#if 0 +"\n"; + s << "xaxis: " << sg.xaxis() << "\n"; + s << "yaxis: " << sg.yaxis() << "\n"; + s << "zaxis: " << sg.zaxis() << "\n"; + + std::vector x = sg.xaxis(); + std::vector y = sg.yaxis(); + std::vector z = sg.zaxis(); + + for ( size_t iz=0 ; iz +#include "histogram.h" + +#include "TH1D.h" +#include "TDirectory.h" + + +class aTH1D; + +TH1D* convert( const aTH1D* h ); +aTH1D* convert( const TH1D* h ); + + +class aTH1D : public histogram { + +public: + + aTH1D() { } + + aTH1D( const histogram& h ) : histogram(h) { } + + aTH1D( const std::string& name, const std::string& title, int nbins, double lo, double hi ) : + histogram( name, nbins, lo, hi ), mtitle(title) { } + + aTH1D( const std::string& name, const std::string& title, int nbins, const double* limits ) : + histogram( name, nbins, limits ), mtitle(title) { } + + virtual ~aTH1D() { } + + double GetBinContent( int i ) const { return my[i-1]; } + double GetBinError( int i ) const { return mye[i-1]; } + + void SetBinContent( int i, double d ) { my[i-1]=d; } + void SetBinError( int i, double d ) { mye[i-1]=d; } + + double GetBinLowEdge( int i ) const { return mxlimits[i-1]; } + double GetBinCenter( int i ) const { return mx[i-1]; } + + double GetBinWidth( int i ) const { return mxlimits[i]-mxlimits[i-1]; } + + int GetNbinsX() const { return int(mx.size()); } + + void Fill( double x, double w=1 ) { fill(x, w ); } + + int FindBin( double x ) const { return index(x)+1; } + + void SetDirectory( int =0 ) const { } + + aTH1D* Clone( const std::string& n="" ) { + aTH1D* h = new aTH1D(*this); + h->SetName( n ); + return h; + } + + void Scale( double d ) { operator*=(d); } + + void Add( const aTH1D& h ) { operator+=(h); } + + void Sumw2() const { } + + void SetName( const std::string& n ) { mname=n; } + std::string GetName() const { return mname; } + + std::string GetTitle() const { return mtitle; } + void SetTitle( const std::string& s ) { mtitle=s; } + + void Reset() { + size_t n = my.size(); + my = std::vector(n,0); + mye = std::vector(n,0); + } + + void Write() const { + TH1D* h = convert( this ); + h->Write(); + } + + +protected: + + std::string mtitle; + +}; + + + +// convertors to and from root histograms + +inline aTH1D* convert( TH1D* h ) { + std::vector limits(h->GetNbinsX()+1); + for ( size_t i=0 ; iGetBinLowEdge(i+1); + aTH1D* hist = new aTH1D( h->GetName(), h->GetTitle(), limits.size()-1, &limits[0] ); + for ( size_t i=0 ; isize() ; i++ ) { + hist->y(i) = h->GetBinContent(i+1); + hist->ye(i) = h->GetBinError(i+1); + } + return hist; +} + + +inline TH1D* convert( const aTH1D* h ) { + TH1D* hist = new TH1D( h->name().c_str(), h->GetTitle().c_str(), h->size(), &(h->xlimits()[0]) ); + for ( size_t i=0 ; isize() ; i++ ) { + hist->SetBinContent(i+1, h->y(i) ); + hist->SetBinError(i+1, h->ye(i) ); + } + return hist; +} + + +inline TH1D* convert( const aTH1D& h ) { + TH1D* hist = new TH1D( h.name().c_str(), h.GetTitle().c_str(), h.size(), &(h.xlimits()[0]) ); + for ( size_t i=h.size() ; i-- ; ) { + hist->SetBinContent(i+1, h.y(i) ); + hist->SetBinError(i+1, h.ye(i) ); + } + return hist; +} + + + +inline std::ostream& operator<<( std::ostream& s, const aTH1D& a ) { + return s << *dynamic_cast(&a); +} + + + + +#endif // ATH1D_H + + + + + + + + + + Index: applgrid/trunk/appl_grid/file_index.h =================================================================== --- applgrid/trunk/appl_grid/file_index.h (revision 0) +++ applgrid/trunk/appl_grid/file_index.h (revision 1895) @@ -0,0 +1,185 @@ +/* emacs: this is -*- c++ -*- */ +/** + ** @file file_index.h + ** + ** @author sutt + ** @date Thu 31 Dec 2020 14:42:12 GMT + ** + ** $Id: file_index.h, v0.0 Thu 31 Dec 2020 14:42:12 GMT sutt $ + ** + ** Copyright (C) 2020 sutt (sutt@cern.ch) + ** + **/ + + +#ifndef FILE_INDEX_H +#define FILE_INDEX_H + +#include +#include +#include +#include + +#include "appl_grid/serialise_base.h" + +class file_index { + +public: + + // histogram error exception + class exception : public std::exception { + public: + exception(const std::string& s) { std::cerr << what() << " " << s << std::endl; }; + virtual const char* what() const throw() { return ""; } + }; + + struct entry { + entry( int sz=0, int off=0 ) : size(sz), offset(off) { } + double size; + double offset; + }; + + +public: + + typedef std::map map_t; + typedef std::map rmap_t; + +public: + + file_index() : mname("index"), mrunning(0) { } + + virtual ~file_index() { } + + std::string name() const { return mname; } + + void add( const std::string& s, int bytes, int offset=-1 ) { + mkeys.push_back( s ); + mmap.insert( map_t::value_type( s, entry( bytes, mrunning ) ) ); + mrmap.insert( rmap_t::value_type( mrunning, s ) ); + if ( offset!=-1 && mrunning!=offset ) std::cerr << "index:: offset mismatch" << std::endl; + mrunning += bytes; + } + + size_t size() const { return mkeys.size(); } + + std::vector::const_iterator begin() const { return mkeys.begin(); } + std::vector::const_iterator end() const { return mkeys.end(); } + + entry find( const std::string& s ) const { + map_t::const_iterator itr = mmap.find( s ); + if ( itr!=mmap.end() ) return itr->second; + return entry(); + } + + + void clear() { + mkeys.clear(); + mmap.clear(); + mrmap.clear(); + } + + + /// serialisation and deserialisation ... + + std::vector serialise() const { + + std::vector s; + + s.push_back( SB::WRITEGUARD ); + + s.push_back(0); // overall size + + SB::serialise( s, name() ); + + s.push_back( mkeys.size() ); + + for ( size_t i=0 ; i& v ) { + + mname.clear(); + mkeys.clear(); + mmap.clear(); + mrmap.clear(); + + std::vector::const_iterator itr = v.begin(); + + if ( (*itr++) != SB::WRITEGUARD ) throw exception("file_index: read error"); + + // int ntotal = (*itr++); + (*itr++); + + SB::deserialise( itr, mname ); + + size_t nkeys = (*itr++); + + for ( size_t i=0 ; i mkeys; + + map_t mmap; + rmap_t mrmap; + + int mrunning; + +}; + +inline std::ostream& operator<<( std::ostream& s, const file_index::entry& e ) { + return s << "index: 0x" << std::hex << int(e.offset) << std::dec << "\tsize: " << int(e.size); +} + + +inline std::ostream& operator<<( std::ostream& s, const file_index& f ) { + s << "file index: " << f.size() << " keys\n"; + for ( std::vector::const_iterator itr=f.begin() ; itr!=f.end() ; itr++ ) { + std::string space = ""; + if ( itr->size()<31 ) space += "\t"; + if ( itr->size()<7 ) space += "\t"; + s << " Key: " << *itr << space << " " << f.find(*itr) << "\n"; + } + return s; +} + + +#endif // FILE_INDEX_H + + + + + + + + + + Index: applgrid/trunk/appl_grid/binary_search.h =================================================================== --- applgrid/trunk/appl_grid/binary_search.h (revision 0) +++ applgrid/trunk/appl_grid/binary_search.h (revision 1895) @@ -0,0 +1,63 @@ +/** emacs: this is -*- c++ -*- **/ +/** + ** @file binary.h + ** + ** fast binary search algorithm + ** + ** @author sutt + ** @date Tue 29 Dec 2020 15:53:45 GMT + ** + ** $Id: binary.h, v0.0 Tue 29 Dec 2020 15:53:45 GMT sutt $ + ** + ** Copyright (C) 2020 sutt (sutt@cern.ch) + ** + **/ + + +#ifndef BINARY_H +#define BINARY_H + +#include + +inline int find_index( double x, const std::vector& v ) { + + if ( x=v.back() ) return -1; + + int a = 0; + int b = v.size()-1; + + while ( b>(a+1) ) { + + /// fast integer division by two + int c = ( a+b ) >> 1; + + if ( x::const_iterator find_lowerbound( double x, const std::vector& v ) { +// return std::lower_bound::const_iterator>( v.begin(), v.end(), x ); +// } + + +#endif // BINARY_H + + + + + + + + + + Index: applgrid/trunk/appl_grid/serialise_base.h =================================================================== --- applgrid/trunk/appl_grid/serialise_base.h (revision 0) +++ applgrid/trunk/appl_grid/serialise_base.h (revision 1895) @@ -0,0 +1,118 @@ +/** emacs: this is -*- c++ -*- **/ +/** + ** @file serialise_base.h + ** + ** @author sutt + ** @date Tue 29 Dec 2020 13:57:14 GMT + ** + ** $Id: serialise_base.h, v0.0 Tue 29 Dec 2020 13:57:14 GMT sutt $ + ** + ** Copyright (C) 2020 sutt (sutt@cern.ch) + ** + **/ + +#ifndef SERIALISE_BASE_H +#define SERIALISE_BASE_H + +#include +#include +#include + +#include + + +namespace SB { + + /// serialisation + + const double WRITEGUARD = 1234567890123456; + + inline void serialise( std::vector& strm, const double& d ) { + strm.push_back( d ); + } + +#if 0 + /// for some reason it cannot resolve the recursive template + /// for string and double + template + inline void serialise( std::vector >& strm, const std::vector& v ) { + strm.push_back( v.size() ); + for ( size_t i=0 ; i& strm, const std::string& s ) { + strm.push_back( s.size() ); + for ( size_t i=0 ; i& strm, const std::vector& v ) { + strm.push_back( v.size() ); + strm.insert( strm.end(), v.begin(), v.end() ); + } + + inline void serialise( std::vector& strm, const std::vector& v ) { + strm.push_back( v.size() ); + for ( size_t i=0 ; i > for some reason + template + inline void serialise( std::vector& strm, const std::vector >& v ) { + strm.push_back( v.size() ); + for ( size_t i=0 ; i::const_iterator& itr, double& v ) { + v = (*itr++); + } + + inline void deserialise( std::vector::const_iterator& itr, std::vector& v ) { + size_t n = (*itr++); + v.clear(); + v.resize(n); + for ( size_t i=0 ; i::const_iterator& itr, std::vector >& v ) { + size_t n = (*itr++); + v.clear(); + v.resize(n); + for ( size_t i=0 ; i::const_iterator& itr, std::string& s ) { + s.clear(); + size_t n = (*itr++); + for ( size_t i=0 ; i::const_iterator& itr, std::vector& v ) { + size_t n = (*itr++); + v.clear(); + v.resize(n); + for ( size_t i=0 ; i::const_iterator& itr, std::vector >& v ) { + size_t n = (*itr++); + v.clear(); + v.resize(n); + for ( size_t i=0 ; i + +#include "serialise_base.h" + +template +class stream_vector { + +public: + + /// constructors ... + + stream_vector() { } + + stream_vector( const std::string& name ) : mname(name) { } + + stream_vector( const std::string& name, const std::vector& v ) : mname(name), mpayload(v) { } + + /// constructer from a stream + stream_vector( const std::vector& v ) { deserialise( v ); } + + stream_vector( const stream_vector& sv ) : + mname(sv.mname), + mpayload(sv.mpayload) + { } + + stream_vector& operator=( stream_vector&& sv ) { + std::swap( mname, sv.mname ); + std::swap( mpayload, sv.mpayload ); + return *this; + } + + stream_vector& operator=( const stream_vector& sv ) { + mname = sv.mname; + mpayload = sv.mpayload; + return *this; + } + + + /// destructor + + virtual ~stream_vector() { } + + + + /// accessors ... + + std::string name() const { return mname; } + + size_t size() const { return mpayload.size(); } + + std::vector& payload() { return mpayload; } + const std::vector& payload() const { return mpayload; } + + T operator[]( int i ) const { return mpayload.at(i); } + T& operator[]( int i ) { return mpayload.at(i); } + + T at( int i ) const { return mpayload.at(i); } + T& at( int i ) { return mpayload.at(i); } + + void push_back( const T& t ) { mpayload.push_back(t); } + void add( const T& t ) { mpayload.push_back(t); } + + typename std::vector::iterator begin() { return mpayload.begin(); } + typename std::vector::iterator end() { return mpayload.end(); } + + typename std::vector::const_iterator begin() const { return mpayload.begin(); } + typename std::vector::const_iterator end() const { return mpayload.end(); } + + + /// serialise ... + + std::vector serialise() const { + + std::vector s; + + s.reserve( 5 + name().size() + payload().size() ); + + /// 5, ie header + trailer + overall size + name size, + /// and number of entries in gthe payload + + s.push_back( SB::WRITEGUARD ); + s.push_back(0); // overall size - updated later + + SB::serialise( s, name() ); + SB::serialise( s, payload() ); + + s.push_back( SB::WRITEGUARD ); + + s[1] = s.size(); + + return s; + } + + + /// deserialise ... + + void deserialise( const std::vector& v ) { + + std::vector::const_iterator itr = v.begin(); + + if ( (*itr++) != SB::WRITEGUARD ) std::cerr << "error reading vector" << std::endl; + + + /// don't keep track of the total object size here + /// just yet ... + // size_t ntotal = (*itr++); + (*itr++); + + SB::deserialise( itr, mname ); + + SB::deserialise( itr, mpayload ); + + double var = (*itr++); + + // if ( ntotal != (5+mname.size()+mpayload.size()) ) std::cerr << "count mismatch: " << name() << std::endl; + + if ( var != SB::WRITEGUARD ) std::cerr << "error reading vector " << name() << std::endl; + } + + +private: + + std::string mname; + + std::vector mpayload; + +}; + +#include "appl_grid/vector_stream.h" + + +template +std::ostream& operator<<( std::ostream& s, const stream_vector& sv ) { + s << "key: " << sv.name() << std::endl; + for ( size_t i=0 ; i + +template +std::ostream& operator<<(std::ostream& s, const std::vector& v) { + for ( unsigned i=0 ; i