Gaudi Framework, version v23r2

Home   Generated: Thu Jun 28 2012

RootTool.h

Go to the documentation of this file.
00001 
00002 /*
00003 * Gaudi namespace declaration
00004 */
00005 namespace Gaudi {
00006 
00016   class RootTool : virtual public RootDataConnection::Tool {
00017   public:
00019     RootTool(RootDataConnection* con) { c = con; }
00021     virtual TBranch* getBranch(CSTR section, CSTR branch_name) {
00022       std::string n = branch_name+".";
00023       for(int i=0, m=n.length()-1; i<m; ++i) if ( !isalnum(n[i]) ) n[i]='_';
00024       TTree* t = c->getSection(section);
00025       TBranch* b = t ? t->GetBranch(n.c_str()) : 0;
00026       if ( !b ) b = t ? t->GetBranch(branch_name.c_str()) : 0;
00027       if (  b ) b->SetAutoDelete(kFALSE);
00028       return b;
00029     }
00031     virtual int loadRefs(CSTR section, CSTR cnt, unsigned long entry, RootObjectRefs& refs)   {
00032       TBranch* b = getBranch(section,cnt+"#R");
00033       RootObjectRefs* prefs = &refs;
00034       if ( b ) {
00035         b->SetAddress(&prefs);
00036         int nb = b->GetEntry(entry);
00037         if ( nb >= 1 )   {
00038           const MergeSections& ms = c->mergeSections();
00039           if ( !ms.empty() ) {
00040             MsgStream& msg = msgSvc();
00041             msgSvc() << MSG::VERBOSE;
00042             pair<const RootRef*,const ContainerSection*> ls = c->getMergeSection(cnt,entry);
00043             if ( ls.first )  {
00044               if ( ls.first->dbase >= 0 ) {
00045                 // Now patch the references and links 'en block' to be efficient
00046                 // First the leafs from the TES
00047                 if ( msg.isActive() ) {
00048                   msg << "Refs: LS [" << entry << "] -> " 
00049                     << ls.first->dbase << "," << ls.first->container 
00050                     << "," << ls.first->link 
00051                     << "," << ls.first->entry
00052                     << endmsg;
00053                 }
00054                 for(size_t j=0, n=refs.refs.size(); j<n; ++j)  {
00055                   RootRef& r = refs.refs[j];
00056                   if ( r.entry>= 0 && r.dbase >= 0 ) {
00057                     int db = r.dbase + ls.first->dbase;
00058                     if ( c->getDb(db) == c->fid() ) {
00059                       if ( r.dbase      ) r.dbase     += ls.first->dbase;
00060                       if ( r.container  ) r.container += ls.first->container;
00061                       if ( r.link       ) r.link      += ls.first->link;
00062                       const string& rc = c->getCont(r.container);
00063                       MergeSections::const_iterator k = ms.find(rc);
00064                       if ( k != ms.end() )   {
00065                         const ContainerSections& cs = (*k).second;
00066                         r.entry = ( ls.first->entry >= 0 && ls.first->entry < (int)cs.size() )
00067                           ? cs[ls.first->entry].start + r.entry : -1;
00068                         if ( msg.isActive() ) {
00069                           msg << "Add link [" << r.entry 
00070                             << "," << ls.first->entry 
00071                             << "," << ls.first->container
00072                             << "," << r.container
00073                             << "," << r.entry
00074                             << "] to -> " << rc << endmsg;
00075                         }
00076                       }
00077                       else {
00078                         msg << MSG::WARNING << c->fid() << "  [" << c->getDb(db) 
00079                           << "] Evt:" << entry << " Invalid link to " << rc << endmsg;
00080                         msg << MSG::VERBOSE;
00081                       }
00082                     }
00083                   }
00084                 }
00086                 for(vector<int>::iterator i=refs.links.begin(); i!=refs.links.end();++i) {
00087                   (*i) += ls.first->link;
00088                 }
00089               }
00090               return nb;
00091             }
00092             return -1;
00093           }
00094           return nb;
00095         }
00096       }
00097       return -1;
00098     }
00100     void addParam(ParamMap& c, char* p) {
00101       char* q = strchr(p,'=');
00102       if ( q ) {
00103         *q = 0;
00104         c.push_back(make_pair(p,++q));
00105       }
00106     }
00108     void addEntry(StringVec& c, char* val) {      c.push_back(val);    }
00110     template <class C, class F> StatusCode readBranch(TTree* t, const char* nam, C& v, F pmf) {
00111       char text[2048];
00112       TBranch* b = t->GetBranch(nam);
00113       if ( b ) {
00114         TLeaf* l = b->GetLeaf(nam);
00115         if ( l ) {
00116           b->SetAddress(text);
00117           msgSvc() << MSG::VERBOSE;
00118           for(Long64_t i=0, n=b->GetEntries(); i<n; ++i) {
00119             if ( b->GetEntry(i)>0 ) {
00120               char* p = (char*)l->GetValuePointer();
00121               msgSvc() << "Add Value[" << b->GetName() << "]:" << p << endmsg;
00122               (this->*pmf)(v,p);
00123             }
00124           }
00125           return StatusCode::SUCCESS;
00126         }
00127       }
00128       msgSvc() << MSG::ERROR << "Failed to read '" << nam << "' table." << endmsg;
00129       return StatusCode::FAILURE;
00130     }
00132     bool get(const string& dsc, pair<string,ContainerSection>& e) {
00133       if ( dsc != "[END-OF-SECTION]" ) { 
00134         size_t id1 = dsc.find("[CNT=");
00135         size_t id2 = dsc.find("[START=");
00136         size_t id3 = dsc.find("[LEN=");
00137         if ( id1 != string::npos && id2 != string::npos && id3 != string::npos ) {
00138           string tmp;
00139           string cnt = dsc.substr(id1+5, id2-1-5);
00140           int section_start  = ::atoi((tmp=dsc.substr(id2+7,id3-id2-8)).c_str());
00141           int section_length = ::atoi((tmp=dsc.substr(id3+5,dsc.find("]",id3+5)-id3-5)).c_str());
00142           e.first = cnt;
00143           e.second = ContainerSection(section_start,section_length);
00144           return true;
00145         }
00146       }
00147       e.first = "";
00148       e.second = ContainerSection(-1,-1);
00149       return false;
00150     }
00152     void analyzeMergeMap(StringVec& tmp) {
00153       StringVec::const_iterator i;
00154       LinkSections&  ls = linkSections();
00155       MergeSections& ms = mergeSections();
00156       pair<string,ContainerSection> e;
00157       MsgStream& msg = msgSvc();
00158       RootRef r;
00159       int cnt = 0;
00160       ls.clear();
00161       ms.clear();
00162       msg << MSG::VERBOSE;
00163       r.dbase = r.container = r.link = r.clid = r.svc = r.entry = 0;
00164       for(i=tmp.begin(); i!=tmp.end();++i) {
00165         if ( get(*i,e) ) {
00166           msg << "Added Merge Section:" << e.first << endmsg;
00167           ms[e.first].push_back(e.second);
00168           if (      e.first == "Links" )
00169             r.link      = e.second.start;
00170           else if ( e.first == "Containers" )
00171             r.container = e.second.start;
00172           else if ( e.first == "Databases" )
00173             r.dbase     = e.second.start;
00174           else if ( e.first == "Params" )
00175             r.svc       = e.second.start;
00176         }
00177         else if ( (*i) == "[END-OF-SECTION]" ) {
00178           r.entry = cnt;
00179           if ( msg.isActive() ) {
00180             msg << "Link Section [" << r.entry << "," << ls.size()
00181               << "] -> D:" << r.dbase
00182               << " C:" << r.container
00183               << " L:" << r.link 
00184               << " P:" << r.svc
00185               << endmsg;
00186           }
00187           ls.push_back(r);
00188           cnt++;
00189         }
00190       }
00191     }
00193     StatusCode readRefs()  {
00194       TTree* t = (TTree*)c->file()->Get("Sections");
00195       StringVec tmp;
00196       if ( t && !readBranch(t,  "Sections",  tmp,       &RootTool::addEntry).isSuccess() )
00197         return StatusCode::FAILURE;
00198       else if ( refs() ) {
00199         analyzeMergeMap(tmp);
00200         if ( !readBranch(refs(),"Databases", dbs(),     &RootTool::addEntry).isSuccess() )
00201           return StatusCode::FAILURE;
00202         if ( !readBranch(refs(),"Containers",conts(),   &RootTool::addEntry).isSuccess() )
00203           return StatusCode::FAILURE;
00204         if ( !readBranch(refs(),"Links",     links(),   &RootTool::addEntry).isSuccess() )
00205           return StatusCode::FAILURE;
00206         if ( !readBranch(refs(),"Params",    params(),  &RootTool::addParam).isSuccess() )
00207           return StatusCode::FAILURE;
00208         return StatusCode::SUCCESS;
00209       }
00210       return StatusCode::FAILURE;
00211     }
00213     string getEntry(const string& c)              {  return c;                    }
00215     string getParam(const pair<string,string>& p) {  return p.first+"="+p.second; }
00217     template <class C, class F> StatusCode saveBranch(const char* nam, C& v, F pmf) {
00218       Long64_t i, n;
00219       string val, typ = nam;
00220       StatusCode sc = StatusCode::SUCCESS;
00221       TDirectory::TContext ctxt(c->file());
00222       TBranch* b = refs()->GetBranch(nam);
00223       if ( !b ) b = refs()->Branch(nam,0,(typ+"/C").c_str());
00224       if ( b ) {
00225         for(i=b->GetEntries(), n=Long64_t(v.size()); i<n; ++i) {
00226           val = (this->*pmf)(v[size_t(i)]);
00227           b->SetAddress((char*)val.c_str());
00228           msgSvc() << MSG::VERBOSE << "Save Value[" << b->GetName() << "]:" << val << endmsg;
00229           if ( b->Fill() <= 1) sc = StatusCode::FAILURE;
00230         }
00231         return sc;
00232       }
00233       return StatusCode::FAILURE;
00234     }
00236     StatusCode saveRefs() {
00237       if ( refs() ) {
00238         if ( !saveBranch("Databases", dbs(),   &RootTool::getEntry).isSuccess() )
00239           return StatusCode::FAILURE;
00240         if ( !saveBranch("Containers",conts(), &RootTool::getEntry).isSuccess() )
00241           return StatusCode::FAILURE;
00242         if ( !saveBranch("Links",     links(), &RootTool::getEntry).isSuccess() )
00243           return StatusCode::FAILURE;
00244         if ( !saveBranch("Params",    params(),&RootTool::getParam).isSuccess() )
00245           return StatusCode::FAILURE;
00246         return StatusCode::SUCCESS;
00247       }
00248       return StatusCode::FAILURE;
00249     }
00250   };
00251 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Defines

Generated at Thu Jun 28 2012 23:27:30 for Gaudi Framework, version v23r2 by Doxygen version 1.7.2 written by Dimitri van Heesch, © 1997-2004