ROOT logo
//*-- Author : Jochen Markert 28.03.2011

//_HADES_CLASS_DESCRIPTION
////////////////////////////////////////////////////////////////////////////
//
// HLoop is a helper class to allow for fast looping of HADES dsts.
// The categories are maped directly from the Tree and allow partical
// reading to speed up. If Hades object exists, the current event structure
// will be replaced by the one defined by HLoop. If Hades not exists, it
// can be created by HLoop constructor. The Hades eventstructure is important
// if one wants to access the data via gHades->getCurrentEvent()->getCategory(cattype)
// or implicit by Classes like HParticleTrackSorter. See the example below for
// further explanation.
//
//
//
//    HEventHeader* getEventHeader()  // retrieve the pointer to the HADES HEventHeader
//    TTree*        getTree       ()  // retrieve the pointer to the HADES tree
//    TChain*       getChain      ()  // retrieve TChain pointer
//    Long64_t      getEntries    ()  // get Number of events in the chain
//
//--------------------------------------------------------
// example:
//
// #include "hloop.h"
// #include "hcategory.h"
// #include "hzip.h"
// #include "heventheader.h"
// #include "hgeantmedia.h"
// #include "hparticlecand.h"
// #include "hparticletracksorter.h"
//
// #include "TIterator.h"
//
// #include <iostream>
// using namespace std;
// void myLoop()
// {
//
//    //---------------LOOP CONFIG----------------------------------------------------------------
//    Bool_t createHades = kTRUE;  // kTRUE = create HADES new
//    HLoop* loop = new HLoop(createHades); // create HADES (needed if one wants to use gHades)
//    loop.setTreeCacheSize(8000000);//8MB (default is 8MB, values <= 0 will disable the cache)
//    //--------------------------------------------------------------------------------------------
//    loop->readSectorFileList("selection.list",kFALSE,kTRUE);  // hldname_without_evtbuilder 1 1 1 1 1 1
//                                                              // a map of all files to active sector is filled
//                                                              // see functions HLoop::getSectors(Int_t*), Bool_t HLoop::goodSector(Int_t sector),HSectorSelector& getSectorSelector()
//                                                              // for file format description see HSectorSelector
//                                                              // parmeters: default value for files not inside the list (kTRUE all sectors on, kFALSE all sectors off)
//                                                              //            print the filelist during read
//                                                              // If no list is read, all sectors are on by default
//                                                              // reading the files list is optional
//                                                              // lines which uses the feature are marked ##SECTOR_SELECTION## in the following example
//    //--------------------------------------------------------------------------------------------
//
//    // add files : Input sources may be combined
//    loop->addFile ("myFile1.root");
//    loop->addFile ("myFile2.root");  // ....
//    loop->addFilesList("list.txt");  // add all files in ascii file list. list contains 1 root file per line (including path)
//    loop->addFiles("myFiles*.root"); // add all files matching this expression
//    loop->addMultFiles("file1.root,file2.root,file3.root");  add files as comma separated list
//    HZip::makeChain("all_files.zip",loop->getChain()); // add all root files contained in zip file (this file has to be produced by hzip)
//
//    if(!loop->setInput("-*,+HParticleCand")) // global disable "-*" has to be first in the list
//    {                                        // read only one category from file (HEventHeader is allways on)
//        exit(1);
//    }
//    loop->printCategories(); // print status from all categories in event
//    //------------------------------------------------------------------------------------------
//
//    TIterator* iterCand = 0;
//    if (loop->getCategory("HParticleCand")) iterCand = loop->getCategory("HParticleCand")->MakeIterator();
//    HEventHeader* header = loop->getEventHeader();
//    HGeantMedia* media = loop->getGeantMedia();
//    if(media){  // NULL for real data or if not in sim input
//       media->print();
//    }
//
//    //--------------------------CONFIGURATION---------------------------------------------------
//    //At begin of the program (outside the event loop)
//
//    HParticleTrackSorter sorter;
//    //sorter.setDebug();                                            // for debug
//    //sorter.setPrintLevel(3);                                      // max prints
//    //sorter.setRICHMatching(HParticleTrackSorter::kUseRKRICHWindow,4.); // select matching RICH-MDC for selectLeptons() function
//    //sorter.setUseYMatching(kTRUE);                                // use meta cell build in select functions
//    //sorter.setMetaBoundary(3.);                                   // match in 3 mm to meta cell
//    //sorter.setUseBeta(kTRUE);                                     // require beta >0 in build in select functions
//    //sorter.setBetaLeptonCut(0.9);                                 // lower beta cut for lepton select
//    //sorter.setUseFakeRejection(kTRUE);                            // reject fakes in  build in select functions (default kTRUE)
//    //sorter.setIgnoreInnerMDC();                                   // do not reject Double_t inner MDC hits
//    //sorter.setIgnoreOuterMDC();                                   // do not reject Double_t outer MDC hits
//    //sorter.setIgnoreMETA();                                       // do not reject Double_t META hits
//    //sorter.setIgnorePreviousIndex();                              // do not reject indices from previous selctions
//    sorter.init();                                                  // get catgegory pointers etc...
//    //--------------------------------------------------------------------------------------------
//    Int_t nFile = 0;
//
//    //--------------------------------------------------------------------------------------------
//    ##SECTOR_SELECTION##
//    Int_t sectors[6];  // array of active sector  (1 active, 0 off) will be retrive inside loop for each event
//    //--------------------------------------------------------------------------------------------
//
//    for (Int_t i = 0; i < 10; i++) {
//        if(loop->nextEvent(i) <= 0) break;  // get next event. categories will be cleared before
//                                            // if 0 (entry not found) or -1 (Io error) is returned stop the loop
//        //--------------------------------------------------------------------------------------------
//        ##SECTOR_SELECTION##
//        loop->getSectors(sectors);  // get active sectors for current event (somethime arrays are useful)
//        //--------------------------------------------------------------------------------------------
//
//        cout<<i<<"----------------------------------- "<<endl;
//        cout<<"sequence Nr = "<<header->getEventSeqNumber()<<endl; // retrieve full header infos
//        TString filename;
//        if(loop->isNewFile(filename)){ // new file opened from chain ?
//            cout<<"new File found "<<filename.Data()<<endl;
//            nFile++;
//	  }
//
//        sorter.cleanUp();
//        sorter.resetFlags(kTRUE,kTRUE,kTRUE,kTRUE);     // reset all flags for flags (0-28) ,reject,used,lepton
//        Int_t nCandLep     = sorter.fill(HParticleTrackSorter::selectLeptons);   // fill only good leptons
//        Int_t nCandLepBest = sorter.selectBest(HParticleTrackSorter::kIsBestRKRKMETA,HParticleTrackSorter::kIsLepton);
//        Int_t nCandHad     = sorter.fill(HParticleTrackSorter::selectHadrons);   // fill only good hadrons (already marked good leptons will be skipped)
//        Int_t nCandHadBest = sorter.selectBest(HParticleTrackSorter::kIsBestRKRKMETA,HParticleTrackSorter::kIsHadron);
//
//        if(iterCand){
//           iterCand->Reset();
//
//           HParticleCand* cand;
//           while ( (cand = (HParticleCand*)iterCand->Next()) != 0 ){
//              // do some work  ....
//              //--------------------------------------------------------------------------------------------
//              ##SECTOR_SELECTION##
//              if(!loop->goodSector(cand->getSector())) continue;  // skipp all candidates from inactive sectors
//              //--------------------------------------------------------------------------------------------
//          }
//       }
//    }
//    sorter.finalize(); // clean up stuff
//
//   if(createHades) delete gHades; // delete Hades or call destructor of HLoop.
// }
//
//----------------------------------------------------------------
//
// A second example running a task with hloop.
// The task can use parameter containers if HSpectrometer
// and parameterio are specified. The task will be executed
// inside nextEvent() function.
//
// void hloopDst()
//{
//    Bool_t createHades=kTRUE;
//    HLoop* loop = new HLoop(createHades);  // kTRUE : create Hades  (needed to work with standard eventstructure)
//    loop.setTreeCacheSize(8000000);//8MB (default is 8MB, values <= 0 will disable the cache)
//
//    TString asciiParFile = "/misc/kempter/projects/plot_dedx/param_31oct11.txt";
//    TString rootParFile = "";
//
//    Int_t mdcMods[6][4]={
//	{0,0,0,0},
//	{0,0,0,0},
//	{1,1,1,1},
//	{0,0,0,0},
//	{1,1,1,1},
//	{1,1,1,1} };
//
//    HDst::setupSpectrometer    ("aug11_3sec",mdcMods,"rich,mdc,tof,rpc,shower,wall,tbox,start");
//    HDst::setupParameterSources("ascii,oracle",asciiParFile,rootParFile,"now");
//
//    loop->addFile("be1122719211201.hld_dst_aug11.root");
//
//    if(!loop->setInput("")) { exit(1); } // read all categories
//
//    HEventHeader* header = loop->getEventHeader();
//    HGeantMedia* media   = loop->getGeantMedia();    // NULL if not in input file ()
//    if(media){ // NULL if not in input file
//        media->print();
//    }
//
//
//    loop->printCategories();
//    loop->printChain();
//
//    HTaskSet *masterTaskSet = gHades->getTaskSet("real");
//    masterTaskSet->add(new HMdcPlotDeDx()); // some task derived from HReconstructor
//
//    HCategory* catCand = HCategoryManager::getCategory(catParticleCand);
//    if(!catCand) { exit(1); }
//
//    Int_t entries=loop->getEntries();
//    for (Int_t i=0; i < entries; i++) {
//	Int_t nbytes =  loop->nextEvent(i);            // get next event. categories will be cleared before
//	if(nbytes <= 0) { cout<<nbytes<<endl; break; } // last event reached
//
//	if(catCand){
//	    HParticleCand* cand=0;
//	    for(Int_t j=0; j<catCand->getEntries();j++){
//		cand = HCategoryManager::getObject(cand,catCand,j);
//                // do something  ....
//	    } // end loop cand
//	} //catCand
//
//    }// end event loop
//    if(gHades)gHades->finalizeTasks();
//    if(createHades) delete gHades; // delete Hades or call destructor of HLoop.
//}
//
////////////////////////////////////////////////////////////////////////////


#include "hloop.h"

#include "hades.h"
#include "htree.h"
#include "hrecevent.h"
#include "hpartialevent.h"
#include "heventheader.h"
#include "hgeantheader.h"
#include "hcategory.h"
#include "htool.h"
#include "htime.h"
#include "hrootsource.h"
#include "hruntimedb.h"
#include "htaskset.h"

#include "TDirectory.h"
#include "TROOT.h"
#include "TGlobal.h"
#include "TList.h"
#include "TKey.h"
#include "TObject.h"
#include "TString.h"
#include "TObjString.h"
#include "TObjArray.h"
#include "TCollection.h"
#include "TOrdCollection.h"
#include "TChainElement.h"
#include "TSystem.h"



#include <iostream>
#include <fstream>
#include <iomanip>

ClassImp(HLoop)
HLoop *gLoop=0;

HLoop::HLoop(Bool_t createHades)
{
    // if createHades == kTRUE (default=kFALSE) Hades will be created and event structure set
    //                == kFALSE Hades will not be created, but if exists
    //                   the event structure will be replaced by the one defined in HLoop
    // if Hades does not exists and is not created new, the event structure will be
    // available only via HLoop and not via gHades
    gLoop = this;
    fChain        = new TChain("T");
    fHead         = NULL;
    fTree         = NULL;
    fFileCurrent  = NULL;
    fGeantMedia   = NULL;
    fMaxEntries   = 0;
    fCurrentEntry = 0;
    fCurrentName  ="unset";
    fIsNewFile    = kTRUE;
    fIsSkipped    = kFALSE;
    fUseTaskSet   = kFALSE;
    fHasCreatedHades = kFALSE;
    fRefID        = -1;
    fFirstEvent   = kTRUE;

    fTreeCacheSize        = 8000000; //8Mb
    fTreeCacheDefaultSize = 8000000; //8Mb;
    fIsCacheSet           = kFALSE;

    fRecEvent     = new HRecEvent();
    for(Int_t s=0;s<6;s++) fsectors[s]=1;

    if(createHades){
	if(gHades == 0) {
            fHasCreatedHades = kTRUE;
	    Hades* hades = new Hades();
            hades->setEvent(fRecEvent);
            Info("HLoop()","HADES created and event structure set !" );
	} else {
	    Error("HLoop()","HADES already exists, but should be created new! Check your macro ..." );
            exit(1);
	}

    } else {
	if(gHades != 0) {
	    gHades->setEvent(fRecEvent);
            Info("HLoop()","HADES already exists, event structure replaced !" );
	} else {
            Info("HLoop()","HADES does not exists! No functions via gHades will be available!" );
	}
    }
}

HLoop::~HLoop()
{
    delete fChain;
    if(fHasCreatedHades && gHades) delete gHades;
    if(fGeantMedia) delete fGeantMedia;
}

void HLoop::setTreeCacheSize(Long64_t cs )
{
    // Configure the read of the Root Tree:
    // default is 8MB to minimize the seek/read on the file system
    // values <= 0 will disable the cache (for comparison)
    // THis function has to be called before setInput().

    fTreeCacheSize = cs;
    fIsCacheSet    = kTRUE;
}


Bool_t HLoop::addFile (TString infile)
{
    // add a single root file to the chain.
    // It will be checked if the file exists.

    TObjArray* elements = fChain->GetListOfFiles();
    Int_t nfiles        = elements->GetEntries();

    if(gSystem->AccessPathName(infile.Data()) == 0) {
	cout<<setw(4)<<nfiles<<" add file : "<<infile.Data()<<endl;
	fChain->AddFile(infile.Data());
        return kTRUE;
    } else {
        Error("addFile()","File = %s not found! Will be skipped .... ",infile.Data());
        return kFALSE;
    }

}

Bool_t HLoop::addFiles(TString expression){

    // add all files matching the expression to the chain

    TObjArray* arr = HTool::glob(expression);
    Int_t n = arr->GetEntries();
    if(n == 0) {
	Warning("addFiles()","No file matching expression '%s' found !",expression.Data());
    } else {
	for(Int_t i = 0; i < n; i ++){
	    TString name = ((TObjString*)arr->At(i))->GetString();
	    addFile(name);
	}
    }
    delete arr;
    return kTRUE;

}
Bool_t HLoop::addFilesList(TString filelist)
{
    // add all files in list filelist. The list should contain
    // 1 root file per line (including path). empty lines and
    // commented lines (starting with "#") will be ignored.
    // Return kFALSE if list files does not exist or and error
    // occured while reading the file

    if(gSystem->AccessPathName(filelist.Data()) == 0 )
    {
	ifstream in;
	in.open(filelist.Data());
	TString name;
        Char_t line[10000];

	while(!in.eof()){
	    in.getline(line,10000);
	    name=line;
	    name.ReplaceAll(" ","");   // remove spaces

	    if(!in.fail()              &&
	       name.CompareTo("")  !=0 &&  // skip empty lines
	       name.BeginsWith("#")==0     // skip commented lines
	      )
	    {
		addFile(name);
	    }
	}
        in.close();

       return kTRUE;
    } else {
       Error("addFilesList()","File = %s not found! Will be skipped .... ",filelist.Data());
       return kFALSE;
    }
}


Bool_t HLoop::addMultFiles(TString commaSeparatedList)
{
    // add multiple files as comma separated file list : "file1,file2,file3"

    TObjArray* arr = commaSeparatedList.Tokenize(",");
    Int_t n = arr->GetEntries();
    if(n == 0) {
	Warning("addMultFiles()","No file in input string '%s' found !",commaSeparatedList.Data());
    } else {
	for(Int_t i = 0; i < n; i ++){
	    TString name = ((TObjString*)arr->At(i))->GetString();
	    if(!addFile(name)){
		delete arr;
		return kFALSE;
	    }
	}
    }
    delete arr;
    return kTRUE;
}


Bool_t HLoop::setInput(TString readCategories)
{
    // read in the given file structure. Categories can be
    // disabled/enabled by the string readCategories. This
    // string should contain a comma separated list. By default
    // all categories are enabled. The format should be
    // like (example 1) "-*,+HMdcCal1,+HParticleCand"
    // (switch off all categories and enable just the to given
    // ones. The global disable -* has to be the first in the
    // list. Or (example 2) "-HMdcRaw,-HRichHit" (enable all
    // categories, but not the give ones).  For simulation
    // or real data one has to use the correct corresponding
    // names (expample : HParticleCandSim / HParticleCand)


    TObjArray* fList  = fChain->GetListOfFiles();

    if(fList->GetEntries() == 0){
        Error("setInput()","No File in List!");
	return kFALSE;
    }



    TChainElement* el = (TChainElement*)fList->At(0);
    TString file = el->GetTitle();

    TFile* fIn = (TFile*)gROOT->GetListOfFiles()->FindObject(file.Data());
    if (!fIn) {
	fIn = new TFile(file.Data());
	if (!fIn) {
	    Error("setFile()","Could not open infile %s",file.Data());
	    return kFALSE;
	}
    }
    fIn->cd();

    fChain->LoadTree(0);

    fTree = (HTree*) fChain->GetTree();
    if(!fTree) {
	Error("setInput()","Could not find Tree in infile %s",file.Data());
	return kFALSE;
    }

    fGeantMedia = (HGeantMedia*) fIn->Get("GeantMedia");


    fHead = fRecEvent->getHeader();


    TBranch* brhead = fTree->GetBranch("EventHeader.");
    if(!brhead) {
	Error("setFile()","Eventheader Branch not found !");
	return kFALSE;
    } else {
	fChain->SetBranchAddress("EventHeader.",&fHead);
	fChain->SetBranchStatus("EventHeader.",1);
	fChain->SetBranchStatus("EventHeader.*",1);
    }

    fChain->SetBranchStatus("Event.",0);
    fChain->SetBranchStatus("Event.*",0);

    TList* dirkeys = gDirectory->GetListOfKeys();

    for(Int_t d = 0; d < dirkeys->GetEntries(); d ++){

	TKey* kdir = (TKey*)dirkeys->At(d);
	TString dirname = kdir->GetName();

	if(dirname.BeginsWith("dir"))
	{

	    TString partialEvent = dirname;
	    partialEvent.ReplaceAll("dir","");

	    if(fIn->cd(dirname.Data()) )
	    {
		TList* catkeys = gDirectory->GetListOfKeys();

		for(Int_t c = 0; c <catkeys->GetEntries(); c ++)
		{
		    fIn->cd(dirname.Data());
		    TKey*   kcat      = (TKey*) catkeys->At(c);
		    TString catname   = kcat->GetName();
		    TString classname = kcat->GetClassName();

		    if(classname == "HPartialEvent") {


			fPartial[catname] = (HPartialEvent*) gDirectory->Get(catname.Data());
                        if(!fPartial[catname]) {
			    Error("setInput()","NULL pointer for partial event %s retrieved!",catname.Data()) ;
			    return kFALSE;
			}

			TBranch* br = fTree->GetBranch(Form("%s.",catname.Data()));
			if(!br) {
			    Error("setInput()","Branch not found : %s !",Form("%s.",catname.Data()));
			    return kFALSE;
			} else {

			    fChain->SetBranchAddress(Form("%s.",catname.Data()),&fPartial[catname]);
			    fChain->SetBranchStatus (Form("%s.",catname.Data()),1);
			    fChain->SetBranchStatus (Form("%s.*",catname.Data()),1);
			}

			fRecEvent->addPartialEvent(fPartial[catname]);

			continue;
		    }

		    //------------------------------------------------------------
		    if(classname == "HLinearCategory" || classname == "HMatrixCategory")
		    {
			fEvent[catname] = (HCategory*) gDirectory->Get(catname.Data());
			if(!fEvent[catname]) {
			    Error("setInput()","NULL pointer for category %s retrieved!",catname.Data()) ;
			    return kFALSE;
			}
			fStatus  [catname] = 0;
                        fPartialN[catname] = partialEvent;

			TBranch* br = fTree->GetBranch(Form("%s.",catname.Data()));
			if(!br) {
			    Error("setInput()","Branch not found : %s !",Form("%s.",catname.Data()));
			    return kFALSE;
			} else {
			    fChain->SetBranchAddress(Form("%s.",catname.Data()),&fEvent[catname]);
			    fChain->SetBranchStatus (Form("%s.",catname.Data()),1);
			    fChain->SetBranchStatus (Form("%s.*",catname.Data()),1);
			    fStatus[catname] = 1;
			}
		    }
		    //------------------------------------------------------------

		} // end cat keys
	    } else {
		Error("setInput()","Could no enter dir %s !",dirname.Data());
		return kFALSE;
	    }
	}
    }
    fIn->cd();

    fMaxEntries = fTree->GetEntries();
    TString readCat = readCategories;
    readCat.ReplaceAll(" ","");

    TObjArray* catNames = readCat.Tokenize(",");
    Int_t n = catNames->GetEntries();

    //------------------------------------------------------------

    for(Int_t i = 0; i < n; i ++){
	TString name = ((TObjString*)catNames->At(i))->GetString();

	//------------------------------------------------------------
	if(i == 0){
	    if(name.CompareTo("-*") == 0){
		setStatus("-*", 0);
		continue;
	    }
	} else {
	    if(name.Contains("*") != 0){
		Warning("setInput()","Wild cards are only supported for '-*' in the first position of the string! Will ignore this one!");
		continue;
	    }
	}
	//------------------------------------------------------------

	//------------------------------------------------------------
	if(name.BeginsWith("+")) { // switch on
	    name.ReplaceAll("+","");
	    if(getCategory(name)) { // exists in input
		setStatus(name, 1);
	    }
	}
	//------------------------------------------------------------

	//------------------------------------------------------------
	if(name.BeginsWith("-")) { // switch off
	    name.ReplaceAll("-","");
	    if(getCategory(name)) { // exists in input
		setStatus(name, 0);
	    }
	}
	//------------------------------------------------------------

    }

    //------------------------------------------------------------
    if(n == 0) {
	// no categories defined in input.
	// by default switch all on
	map<TString,HCategory*>::iterator iter;
	for( iter = fEvent.begin(); iter != fEvent.end(); ++iter ) {
	    addCatName(iter->first,getCategory(iter->first)->getCategory());
	}
    } else {

	map<TString,Int_t>::iterator iter;
	for( iter = fStatus.begin(); iter != fStatus.end(); ++iter ) {
	    if(iter->second == 1 ) {
		addCatName(iter->first,getCategory(iter->first)->getCategory());
	    }
	}
    }
    //------------------------------------------------------------


    fIn->Close();

    printChain();
    Long64_t cs = fChain->GetCacheSize();
    if(!fIsCacheSet ){
	cs = fTreeCacheDefaultSize;
    } else {
        cs = fTreeCacheSize;
    }

    if(fTreeCacheSize > 0) Info("setInput()","Setting Tree Cache size = %lld byte.",cs);
    else                   Info("setInput()","Tree cache is disabled");

    if(cs > 0){
	fChain->SetCacheSize(cs);
	//fChain->AddBranchToCache("*",kTRUE);
	//fChain->StopCacheLearningPhase();
    }
    return kTRUE;
}

void HLoop::printCategories()
{
    // print all categories found in input + status

    cout<<"--------------------------------------"<<endl;
    cout<<"current event :"<<endl;
    map<TString,HCategory*>::iterator iter;
    for( iter = fEvent.begin(); iter != fEvent.end(); ++iter ) {
	cout << setw(20)<<iter->first.Data() << " = " << iter->second <<" status = "<<fStatus[iter->first]<< endl;
    }
    cout<<"--------------------------------------"<<endl;
}

void HLoop::printChain()
{
    // print all files in the chain


    TObjArray* elements = fChain->GetListOfFiles();
    cout<<"--------------------------------------"<<endl;
    cout<<" Analize Chain : .... Can take while  for N file = "<<elements->GetEntries()<<endl;
    fMaxEntries = fChain->GetEntries();
    cout<<"N total entries = "<<fMaxEntries<<endl;
    for(Int_t i = 0; i < elements->GetEntries(); ++i ) {
        TChainElement* element = (TChainElement*)elements->At(i);
	cout <<setw(4)<<i<< setw(80)<<element->GetTitle() << ", N entries = " << element->GetEntries()<< endl;
    }
    cout<<"--------------------------------------"<<endl;
}

void HLoop::printBranchStatus()
{
    // print the branch status of the chain

    cout<<"--------------------------------------"<<endl;
    cout<<"Branch status :"<<endl;

    TObjArray* branches = fChain->GetListOfBranches();
    for(Int_t i = 0; i < branches ->GetEntries(); ++i ) {
        TBranch* br = (TBranch*) branches->At(i);
	cout <<setw(4)<<i<< setw(30)<<br->GetName() << ", status = " <<(Int_t) fChain->GetBranchStatus(br->GetName())<< endl;
    }
    cout<<"--------------------------------------"<<endl;
}

void HLoop::clearCategories()
{
    //  clear all active catgegories of the event

    map<TString,HCategory*>::iterator iter;
    for( iter = fEvent.begin(); iter != fEvent.end(); ++iter ) {
	if(fStatus[iter->first]) iter->second->Clear();
    }
}

HCategory* HLoop::getCategory(TString catname,Bool_t silent)
{
    // return the pointer to given category. NULL if not found

    TString simname  = catname;
    TString realname = catname;

    if(realname.EndsWith("Sim")==1) realname.Replace(realname.Length(),3,"");
    if(simname .EndsWith("Sim")==0) simname+="Sim";


    map<TString,HCategory*>::iterator iter = fEvent.find(realname);
    if( iter == fEvent.end() ){

	iter = fEvent.find(simname);
	if(iter == fEvent.end()){

	    if(!silent)Error("getCategory()","Category \"%s\" not found!",catname.Data());
	    return 0;

	} else  return iter->second;

    } else return iter->second;
}

Bool_t HLoop::getCategoryStatus(TString catname,Bool_t silent)
{
    // return the status of a given category. kFALSE if not found

    map<TString,Int_t>::iterator iter = fStatus.find(catname);
    if( iter == fStatus.end() ){
	if(!silent)Error("getCategory()","Category \"%s\" not found!",catname.Data());
	return kFALSE;
    } else {
	return (Bool_t)fStatus[catname];
    }
}

HPartialEvent* HLoop::getPartialEvent(TString catname,Bool_t silent)
{
    // return the pointer to given partialevent. NULL if not found

    map<TString,HPartialEvent*>::iterator iter = fPartial.find(catname);
    if( iter == fPartial.end() ){
	if(!silent)Error("getPartialEvent()","PartialEvent \"%s\" not found!",catname.Data());
	return 0;
    }
    else return iter->second;
}

HGeantHeader* HLoop::getGeantHeader(Bool_t silent)
{
    // return the pointer to HGeantHader (subheader of PartialEvent Simul).
    // NULL if not found
    HPartialEvent* pEvt = getPartialEvent("Simul",silent);
    HGeantHeader* header = NULL;
    if(pEvt) header = (HGeantHeader*) pEvt->getSubHeader();
    return header;
}

Bool_t HLoop::setStatus(TString catname, Int_t stat )
{
    // set the status of a given category. return kFALSE
    // if something went wrong. if catname="-*" or "+*"
    // all categories of the event are disabled /enabled.

    if(catname.CompareTo("-*") == 0 || catname.CompareTo("+*") == 0 ){
	map<TString,HCategory*>::iterator iter;
	for( iter = fEvent.begin(); iter != fEvent.end(); ++iter ) {
	    if(catname.CompareTo("-*") == 0)  {
		fChain->SetBranchStatus(Form("%s.",iter->first.Data()),0);
		fChain->SetBranchStatus(Form("%s.*",iter->first.Data()),0);
		fStatus[iter->first] = 0;
	    }
	    if(catname.CompareTo("+*") == 0)  {
		fChain->SetBranchStatus(Form("%s.",iter->first.Data()),1);
		fChain->SetBranchStatus(Form("%s.*",iter->first.Data()),1);
		fStatus[iter->first] = 1;
	    }
	}
	return kTRUE;
    }
    if(!getCategory(catname)) {
	return kFALSE;
    }

    TString simname  = catname;
    TString realname = catname;

    if(realname.EndsWith("Sim")==1) realname.Replace(realname.Length(),3,"");
    if(simname .EndsWith("Sim")==0) simname+="Sim";

    map<TString,Int_t>::iterator iter = fStatus.find(realname);

    if( iter == fStatus.end() ){

	iter = fStatus.find(simname);

	if( iter == fStatus.end() ) {
	    Error("getCategory()","Category \"%s\" not found!",catname.Data());
	    return kFALSE;
	} else {

	    fChain->SetBranchStatus(Form("%s.",simname.Data()),stat);
	    fChain->SetBranchStatus(Form("%s.*",simname.Data()),stat);
	    iter->second = stat;
	    return kTRUE;
	}
    }
    else {
	fChain->SetBranchStatus(Form("%s.",realname.Data()),stat);
	fChain->SetBranchStatus(Form("%s.*",realname.Data()),stat);
	iter->second = stat;
	return kTRUE;
    }
    return kFALSE;
}

Int_t HLoop::nextEvent(Int_t iev)
{
    // get next entry form tree. clears all active categories
    // before. Returns number of read bytes
    // (0 = entry not found , -1 = IO error)
    // prints the current file name if a new file
    // is entered in the chain

    if(!fHead) {
	Warning("nextEvent()","File not initialized ! Call setInput(....) before!");
	exit(1);
    }

    fCurrentEntry = iev;

    clearCategories();

    //--------------------------------------------------------
    // read next event from Chain
    Int_t nbytes = fChain->GetEntry(fCurrentEntry);  // 0 = entry not found , -1 = IO error
    if(nbytes == 0) {
	Warning("nextEvent()","Entry %i not found!",(Int_t)fCurrentEntry);
    } else if(nbytes == -1) {
	Warning("nextEvent()","Entry %i read with Io error!",(Int_t)fCurrentEntry);
    }

    //--------------------------------------------------------

    //--------------------------------------------------------
    // check for new file
    TString name = fChain->GetFile()->GetName();
    if(name != fCurrentName) {
	fIsNewFile = kTRUE;
        Int_t entries = 0;
	TObjArray* elements = fChain->GetListOfFiles();
	for(Int_t i = 0; i < elements->GetEntries(); ++i ) {
	    TChainElement* element = (TChainElement*)elements->At(i);
	    if(name.CompareTo(element->GetTitle() ) == 0) {
		entries = element->GetEntries();
                break;
	    }
	}
	fCurrentName = fChain->GetFile()->GetName() ;
	fFileCurrent = fChain->GetFile() ;
	fTree        = fChain->GetTree() ;

	//--------------------------------------------------------
	// sync the sector array from filelist
	// 1. do not try if there has been no filelist read
	// 2. do silent error
	// 3. local fsectors will hold the actual vals
	// 4. fsectors will be filled with default 1 if filename does not match hldname format (2chars+13digits)
        // 5. for correct filenames which are not inside the list fsectors will be filled by the default val set in the readSectorFileList()

	for (Int_t s=0;s<6;s++){ fsectors[s]=1;} // defaults is on
	if(fSectorSelector.getNFiles() > 0){ // no sector file list has been read
	    TString hldname = HTime::stripFileName(fCurrentName,kTRUE,kTRUE); // remove Evtbuild from hldname, do it silent if there is an error
	    if(hldname != "" ) { // filename was compatible
		fSectorSelector.getSectors(fsectors,hldname); // result depends on defaultVal set
	    }
	}
	//--------------------------------------------------------

	cout<<"---------------------------------------------------------"<<endl;
	cout<<name.Data() <<",  sectors : "<<fsectors[0]<<" "<<fsectors[1]<<" "<<fsectors[2]<<" "<<fsectors[3]<<" "<<fsectors[4]<<" "<<fsectors[5]<<" "<<endl;
	cout<<"at entry "<<setw(8)<<fCurrentEntry<<" with entries = "<<entries<<endl;
        cout<<"---------------------------------------------------------"<<endl;

    } else {
       fIsNewFile = kFALSE;
    }
    //--------------------------------------------------------

    if(fFirstEvent){
	if(gHades){
	    //--------------------------------------------------------
	    const TObjArray* l = gHades->getTaskList();

	    HTaskSet* task;
	    TIter iter(l);
	    Int_t ct = 0 ;
	    while((task = (HTaskSet*)iter())) {
		TOrdCollection* list = task->getSetOfTask();
		ct += list->GetEntries();
	    }

	    if(ct > 0) fUseTaskSet = kTRUE;

	    if(fUseTaskSet && gHades->getSetup() == 0 ) {
		Error("","TaskSet used, but Spectrometer not set!"); exit(1);
		fUseTaskSet = kFALSE;
	    }
	    //--------------------------------------------------------

	    if(fUseTaskSet){

		HRootSource* source = new HRootSource();
		source->addFile((Text_t*)fCurrentName.Data());
		source->setCurrentRunId(fHead->getEventRunNumber());
		source->setGlobalRefId(fRefID);
		gHades->setDataSource(source);

		if(!gHades->init(kTRUE)){ Error("","HADES init failed!"); exit(1);}

		HRuntimeDb* rtdb = gHades->getRuntimeDb();

		if (rtdb){
		    if(!rtdb->initContainers(fHead->getEventRunNumber(),
					     fRefID,
					     fCurrentName.Data())
		      ){
			Error("","Rtdb init failed!"); exit(1);
		    }
		}
		gHades->reinitTasks();
		gHades->printDefinedTaskSets ();
                rtdb->disconnectInputs();
	    }
	}
	fFirstEvent = kFALSE;
    }

    //--------------------------------------------------------
    // reinit if new file is opened
    if(fUseTaskSet && !fFirstEvent && gHades && fIsNewFile){
	HRootSource* source = (HRootSource*)gHades->getDataSource();
	if(source){
	    HRuntimeDb* rtdb = gHades->getRuntimeDb();
	    rtdb->reconnectInputs();
	    source->setCurrentRunId(fHead->getEventRunNumber());
	    source->setGlobalRefId(fRefID);

	    if(!rtdb->initContainers(fHead->getEventRunNumber(),
				     fRefID,
				     fCurrentName.Data())
	      ){
		Error("","Rtdb init failed!"); exit(1);
	    }

	    gHades->reinitTasks();
	    rtdb->disconnectInputs();
	}
    }
    //--------------------------------------------------------

    //--------------------------------------------------------
    // run tasks if they are defined
    if(fUseTaskSet && gHades){
	Int_t retVal = gHades->executeTasks();
	if(retVal == kSkipEvent){
	    fIsSkipped = kTRUE;
	    fRecEvent->Clear();
	} else {
            fIsSkipped = kFALSE;
	}
    }
    //--------------------------------------------------------


    return nbytes;
}

Bool_t HLoop::addCatName(TString catname,Short_t catNum)
{
    // add agiven category number to the HADES event structure.
    // this function has to be called after setInput() !!!!

    if(fHead == 0) {
        Error("addCatName()","You have to call setInput() before calling this function !");
	exit(1);
    }

    map<TString,Short_t>::iterator iter = fNameToCat.find(catname);
    if( iter == fNameToCat.end() ){
	if(getCategory(catname,kTRUE)){
            Bool_t ret =  fRecEvent->addCategory((Cat_t)catNum,getCategory(catname),fPartialN[catname].Data());
            fNameToCat[catname] = catNum;
	    return ret;
       } else {
	   Error("addCatName()","Category \"%s\" does not exists in input! Will skip ...",catname.Data());
           return kFALSE;
       }
    } else {
	Error("addCatName()","Category \"%s\" already exists with number = %i!",catname.Data(),fNameToCat[catname]);
        return kFALSE;
    }
    return kTRUE;
}

Bool_t HLoop::isNewFile(TString& name)
{
    // this function can be called after HLoop::nextEvent(Int_t iev)
    // If a new file was opened kTRUE will be returned. The
    // current file name will be returned by name
    name = fCurrentName;
    return fIsNewFile;
}

TObject* HLoop::getFromInputFile(TString name)
{
    // returns an object from root input file
    // return 0 if not existsing, file not opened
    // or anything else strange.
    if(name=="") return 0;

    if(fFileCurrent){
	TDirectory* dold =  gDirectory;

        fFileCurrent->cd();

        TObject* obj = fFileCurrent->Get(name.Data());

	dold->cd();

        return obj;

    } else return 0;
}


 hloop.cc:1
 hloop.cc:2
 hloop.cc:3
 hloop.cc:4
 hloop.cc:5
 hloop.cc:6
 hloop.cc:7
 hloop.cc:8
 hloop.cc:9
 hloop.cc:10
 hloop.cc:11
 hloop.cc:12
 hloop.cc:13
 hloop.cc:14
 hloop.cc:15
 hloop.cc:16
 hloop.cc:17
 hloop.cc:18
 hloop.cc:19
 hloop.cc:20
 hloop.cc:21
 hloop.cc:22
 hloop.cc:23
 hloop.cc:24
 hloop.cc:25
 hloop.cc:26
 hloop.cc:27
 hloop.cc:28
 hloop.cc:29
 hloop.cc:30
 hloop.cc:31
 hloop.cc:32
 hloop.cc:33
 hloop.cc:34
 hloop.cc:35
 hloop.cc:36
 hloop.cc:37
 hloop.cc:38
 hloop.cc:39
 hloop.cc:40
 hloop.cc:41
 hloop.cc:42
 hloop.cc:43
 hloop.cc:44
 hloop.cc:45
 hloop.cc:46
 hloop.cc:47
 hloop.cc:48
 hloop.cc:49
 hloop.cc:50
 hloop.cc:51
 hloop.cc:52
 hloop.cc:53
 hloop.cc:54
 hloop.cc:55
 hloop.cc:56
 hloop.cc:57
 hloop.cc:58
 hloop.cc:59
 hloop.cc:60
 hloop.cc:61
 hloop.cc:62
 hloop.cc:63
 hloop.cc:64
 hloop.cc:65
 hloop.cc:66
 hloop.cc:67
 hloop.cc:68
 hloop.cc:69
 hloop.cc:70
 hloop.cc:71
 hloop.cc:72
 hloop.cc:73
 hloop.cc:74
 hloop.cc:75
 hloop.cc:76
 hloop.cc:77
 hloop.cc:78
 hloop.cc:79
 hloop.cc:80
 hloop.cc:81
 hloop.cc:82
 hloop.cc:83
 hloop.cc:84
 hloop.cc:85
 hloop.cc:86
 hloop.cc:87
 hloop.cc:88
 hloop.cc:89
 hloop.cc:90
 hloop.cc:91
 hloop.cc:92
 hloop.cc:93
 hloop.cc:94
 hloop.cc:95
 hloop.cc:96
 hloop.cc:97
 hloop.cc:98
 hloop.cc:99
 hloop.cc:100
 hloop.cc:101
 hloop.cc:102
 hloop.cc:103
 hloop.cc:104
 hloop.cc:105
 hloop.cc:106
 hloop.cc:107
 hloop.cc:108
 hloop.cc:109
 hloop.cc:110
 hloop.cc:111
 hloop.cc:112
 hloop.cc:113
 hloop.cc:114
 hloop.cc:115
 hloop.cc:116
 hloop.cc:117
 hloop.cc:118
 hloop.cc:119
 hloop.cc:120
 hloop.cc:121
 hloop.cc:122
 hloop.cc:123
 hloop.cc:124
 hloop.cc:125
 hloop.cc:126
 hloop.cc:127
 hloop.cc:128
 hloop.cc:129
 hloop.cc:130
 hloop.cc:131
 hloop.cc:132
 hloop.cc:133
 hloop.cc:134
 hloop.cc:135
 hloop.cc:136
 hloop.cc:137
 hloop.cc:138
 hloop.cc:139
 hloop.cc:140
 hloop.cc:141
 hloop.cc:142
 hloop.cc:143
 hloop.cc:144
 hloop.cc:145
 hloop.cc:146
 hloop.cc:147
 hloop.cc:148
 hloop.cc:149
 hloop.cc:150
 hloop.cc:151
 hloop.cc:152
 hloop.cc:153
 hloop.cc:154
 hloop.cc:155
 hloop.cc:156
 hloop.cc:157
 hloop.cc:158
 hloop.cc:159
 hloop.cc:160
 hloop.cc:161
 hloop.cc:162
 hloop.cc:163
 hloop.cc:164
 hloop.cc:165
 hloop.cc:166
 hloop.cc:167
 hloop.cc:168
 hloop.cc:169
 hloop.cc:170
 hloop.cc:171
 hloop.cc:172
 hloop.cc:173
 hloop.cc:174
 hloop.cc:175
 hloop.cc:176
 hloop.cc:177
 hloop.cc:178
 hloop.cc:179
 hloop.cc:180
 hloop.cc:181
 hloop.cc:182
 hloop.cc:183
 hloop.cc:184
 hloop.cc:185
 hloop.cc:186
 hloop.cc:187
 hloop.cc:188
 hloop.cc:189
 hloop.cc:190
 hloop.cc:191
 hloop.cc:192
 hloop.cc:193
 hloop.cc:194
 hloop.cc:195
 hloop.cc:196
 hloop.cc:197
 hloop.cc:198
 hloop.cc:199
 hloop.cc:200
 hloop.cc:201
 hloop.cc:202
 hloop.cc:203
 hloop.cc:204
 hloop.cc:205
 hloop.cc:206
 hloop.cc:207
 hloop.cc:208
 hloop.cc:209
 hloop.cc:210
 hloop.cc:211
 hloop.cc:212
 hloop.cc:213
 hloop.cc:214
 hloop.cc:215
 hloop.cc:216
 hloop.cc:217
 hloop.cc:218
 hloop.cc:219
 hloop.cc:220
 hloop.cc:221
 hloop.cc:222
 hloop.cc:223
 hloop.cc:224
 hloop.cc:225
 hloop.cc:226
 hloop.cc:227
 hloop.cc:228
 hloop.cc:229
 hloop.cc:230
 hloop.cc:231
 hloop.cc:232
 hloop.cc:233
 hloop.cc:234
 hloop.cc:235
 hloop.cc:236
 hloop.cc:237
 hloop.cc:238
 hloop.cc:239
 hloop.cc:240
 hloop.cc:241
 hloop.cc:242
 hloop.cc:243
 hloop.cc:244
 hloop.cc:245
 hloop.cc:246
 hloop.cc:247
 hloop.cc:248
 hloop.cc:249
 hloop.cc:250
 hloop.cc:251
 hloop.cc:252
 hloop.cc:253
 hloop.cc:254
 hloop.cc:255
 hloop.cc:256
 hloop.cc:257
 hloop.cc:258
 hloop.cc:259
 hloop.cc:260
 hloop.cc:261
 hloop.cc:262
 hloop.cc:263
 hloop.cc:264
 hloop.cc:265
 hloop.cc:266
 hloop.cc:267
 hloop.cc:268
 hloop.cc:269
 hloop.cc:270
 hloop.cc:271
 hloop.cc:272
 hloop.cc:273
 hloop.cc:274
 hloop.cc:275
 hloop.cc:276
 hloop.cc:277
 hloop.cc:278
 hloop.cc:279
 hloop.cc:280
 hloop.cc:281
 hloop.cc:282
 hloop.cc:283
 hloop.cc:284
 hloop.cc:285
 hloop.cc:286
 hloop.cc:287
 hloop.cc:288
 hloop.cc:289
 hloop.cc:290
 hloop.cc:291
 hloop.cc:292
 hloop.cc:293
 hloop.cc:294
 hloop.cc:295
 hloop.cc:296
 hloop.cc:297
 hloop.cc:298
 hloop.cc:299
 hloop.cc:300
 hloop.cc:301
 hloop.cc:302
 hloop.cc:303
 hloop.cc:304
 hloop.cc:305
 hloop.cc:306
 hloop.cc:307
 hloop.cc:308
 hloop.cc:309
 hloop.cc:310
 hloop.cc:311
 hloop.cc:312
 hloop.cc:313
 hloop.cc:314
 hloop.cc:315
 hloop.cc:316
 hloop.cc:317
 hloop.cc:318
 hloop.cc:319
 hloop.cc:320
 hloop.cc:321
 hloop.cc:322
 hloop.cc:323
 hloop.cc:324
 hloop.cc:325
 hloop.cc:326
 hloop.cc:327
 hloop.cc:328
 hloop.cc:329
 hloop.cc:330
 hloop.cc:331
 hloop.cc:332
 hloop.cc:333
 hloop.cc:334
 hloop.cc:335
 hloop.cc:336
 hloop.cc:337
 hloop.cc:338
 hloop.cc:339
 hloop.cc:340
 hloop.cc:341
 hloop.cc:342
 hloop.cc:343
 hloop.cc:344
 hloop.cc:345
 hloop.cc:346
 hloop.cc:347
 hloop.cc:348
 hloop.cc:349
 hloop.cc:350
 hloop.cc:351
 hloop.cc:352
 hloop.cc:353
 hloop.cc:354
 hloop.cc:355
 hloop.cc:356
 hloop.cc:357
 hloop.cc:358
 hloop.cc:359
 hloop.cc:360
 hloop.cc:361
 hloop.cc:362
 hloop.cc:363
 hloop.cc:364
 hloop.cc:365
 hloop.cc:366
 hloop.cc:367
 hloop.cc:368
 hloop.cc:369
 hloop.cc:370
 hloop.cc:371
 hloop.cc:372
 hloop.cc:373
 hloop.cc:374
 hloop.cc:375
 hloop.cc:376
 hloop.cc:377
 hloop.cc:378
 hloop.cc:379
 hloop.cc:380
 hloop.cc:381
 hloop.cc:382
 hloop.cc:383
 hloop.cc:384
 hloop.cc:385
 hloop.cc:386
 hloop.cc:387
 hloop.cc:388
 hloop.cc:389
 hloop.cc:390
 hloop.cc:391
 hloop.cc:392
 hloop.cc:393
 hloop.cc:394
 hloop.cc:395
 hloop.cc:396
 hloop.cc:397
 hloop.cc:398
 hloop.cc:399
 hloop.cc:400
 hloop.cc:401
 hloop.cc:402
 hloop.cc:403
 hloop.cc:404
 hloop.cc:405
 hloop.cc:406
 hloop.cc:407
 hloop.cc:408
 hloop.cc:409
 hloop.cc:410
 hloop.cc:411
 hloop.cc:412
 hloop.cc:413
 hloop.cc:414
 hloop.cc:415
 hloop.cc:416
 hloop.cc:417
 hloop.cc:418
 hloop.cc:419
 hloop.cc:420
 hloop.cc:421
 hloop.cc:422
 hloop.cc:423
 hloop.cc:424
 hloop.cc:425
 hloop.cc:426
 hloop.cc:427
 hloop.cc:428
 hloop.cc:429
 hloop.cc:430
 hloop.cc:431
 hloop.cc:432
 hloop.cc:433
 hloop.cc:434
 hloop.cc:435
 hloop.cc:436
 hloop.cc:437
 hloop.cc:438
 hloop.cc:439
 hloop.cc:440
 hloop.cc:441
 hloop.cc:442
 hloop.cc:443
 hloop.cc:444
 hloop.cc:445
 hloop.cc:446
 hloop.cc:447
 hloop.cc:448
 hloop.cc:449
 hloop.cc:450
 hloop.cc:451
 hloop.cc:452
 hloop.cc:453
 hloop.cc:454
 hloop.cc:455
 hloop.cc:456
 hloop.cc:457
 hloop.cc:458
 hloop.cc:459
 hloop.cc:460
 hloop.cc:461
 hloop.cc:462
 hloop.cc:463
 hloop.cc:464
 hloop.cc:465
 hloop.cc:466
 hloop.cc:467
 hloop.cc:468
 hloop.cc:469
 hloop.cc:470
 hloop.cc:471
 hloop.cc:472
 hloop.cc:473
 hloop.cc:474
 hloop.cc:475
 hloop.cc:476
 hloop.cc:477
 hloop.cc:478
 hloop.cc:479
 hloop.cc:480
 hloop.cc:481
 hloop.cc:482
 hloop.cc:483
 hloop.cc:484
 hloop.cc:485
 hloop.cc:486
 hloop.cc:487
 hloop.cc:488
 hloop.cc:489
 hloop.cc:490
 hloop.cc:491
 hloop.cc:492
 hloop.cc:493
 hloop.cc:494
 hloop.cc:495
 hloop.cc:496
 hloop.cc:497
 hloop.cc:498
 hloop.cc:499
 hloop.cc:500
 hloop.cc:501
 hloop.cc:502
 hloop.cc:503
 hloop.cc:504
 hloop.cc:505
 hloop.cc:506
 hloop.cc:507
 hloop.cc:508
 hloop.cc:509
 hloop.cc:510
 hloop.cc:511
 hloop.cc:512
 hloop.cc:513
 hloop.cc:514
 hloop.cc:515
 hloop.cc:516
 hloop.cc:517
 hloop.cc:518
 hloop.cc:519
 hloop.cc:520
 hloop.cc:521
 hloop.cc:522
 hloop.cc:523
 hloop.cc:524
 hloop.cc:525
 hloop.cc:526
 hloop.cc:527
 hloop.cc:528
 hloop.cc:529
 hloop.cc:530
 hloop.cc:531
 hloop.cc:532
 hloop.cc:533
 hloop.cc:534
 hloop.cc:535
 hloop.cc:536
 hloop.cc:537
 hloop.cc:538
 hloop.cc:539
 hloop.cc:540
 hloop.cc:541
 hloop.cc:542
 hloop.cc:543
 hloop.cc:544
 hloop.cc:545
 hloop.cc:546
 hloop.cc:547
 hloop.cc:548
 hloop.cc:549
 hloop.cc:550
 hloop.cc:551
 hloop.cc:552
 hloop.cc:553
 hloop.cc:554
 hloop.cc:555
 hloop.cc:556
 hloop.cc:557
 hloop.cc:558
 hloop.cc:559
 hloop.cc:560
 hloop.cc:561
 hloop.cc:562
 hloop.cc:563
 hloop.cc:564
 hloop.cc:565
 hloop.cc:566
 hloop.cc:567
 hloop.cc:568
 hloop.cc:569
 hloop.cc:570
 hloop.cc:571
 hloop.cc:572
 hloop.cc:573
 hloop.cc:574
 hloop.cc:575
 hloop.cc:576
 hloop.cc:577
 hloop.cc:578
 hloop.cc:579
 hloop.cc:580
 hloop.cc:581
 hloop.cc:582
 hloop.cc:583
 hloop.cc:584
 hloop.cc:585
 hloop.cc:586
 hloop.cc:587
 hloop.cc:588
 hloop.cc:589
 hloop.cc:590
 hloop.cc:591
 hloop.cc:592
 hloop.cc:593
 hloop.cc:594
 hloop.cc:595
 hloop.cc:596
 hloop.cc:597
 hloop.cc:598
 hloop.cc:599
 hloop.cc:600
 hloop.cc:601
 hloop.cc:602
 hloop.cc:603
 hloop.cc:604
 hloop.cc:605
 hloop.cc:606
 hloop.cc:607
 hloop.cc:608
 hloop.cc:609
 hloop.cc:610
 hloop.cc:611
 hloop.cc:612
 hloop.cc:613
 hloop.cc:614
 hloop.cc:615
 hloop.cc:616
 hloop.cc:617
 hloop.cc:618
 hloop.cc:619
 hloop.cc:620
 hloop.cc:621
 hloop.cc:622
 hloop.cc:623
 hloop.cc:624
 hloop.cc:625
 hloop.cc:626
 hloop.cc:627
 hloop.cc:628
 hloop.cc:629
 hloop.cc:630
 hloop.cc:631
 hloop.cc:632
 hloop.cc:633
 hloop.cc:634
 hloop.cc:635
 hloop.cc:636
 hloop.cc:637
 hloop.cc:638
 hloop.cc:639
 hloop.cc:640
 hloop.cc:641
 hloop.cc:642
 hloop.cc:643
 hloop.cc:644
 hloop.cc:645
 hloop.cc:646
 hloop.cc:647
 hloop.cc:648
 hloop.cc:649
 hloop.cc:650
 hloop.cc:651
 hloop.cc:652
 hloop.cc:653
 hloop.cc:654
 hloop.cc:655
 hloop.cc:656
 hloop.cc:657
 hloop.cc:658
 hloop.cc:659
 hloop.cc:660
 hloop.cc:661
 hloop.cc:662
 hloop.cc:663
 hloop.cc:664
 hloop.cc:665
 hloop.cc:666
 hloop.cc:667
 hloop.cc:668
 hloop.cc:669
 hloop.cc:670
 hloop.cc:671
 hloop.cc:672
 hloop.cc:673
 hloop.cc:674
 hloop.cc:675
 hloop.cc:676
 hloop.cc:677
 hloop.cc:678
 hloop.cc:679
 hloop.cc:680
 hloop.cc:681
 hloop.cc:682
 hloop.cc:683
 hloop.cc:684
 hloop.cc:685
 hloop.cc:686
 hloop.cc:687
 hloop.cc:688
 hloop.cc:689
 hloop.cc:690
 hloop.cc:691
 hloop.cc:692
 hloop.cc:693
 hloop.cc:694
 hloop.cc:695
 hloop.cc:696
 hloop.cc:697
 hloop.cc:698
 hloop.cc:699
 hloop.cc:700
 hloop.cc:701
 hloop.cc:702
 hloop.cc:703
 hloop.cc:704
 hloop.cc:705
 hloop.cc:706
 hloop.cc:707
 hloop.cc:708
 hloop.cc:709
 hloop.cc:710
 hloop.cc:711
 hloop.cc:712
 hloop.cc:713
 hloop.cc:714
 hloop.cc:715
 hloop.cc:716
 hloop.cc:717
 hloop.cc:718
 hloop.cc:719
 hloop.cc:720
 hloop.cc:721
 hloop.cc:722
 hloop.cc:723
 hloop.cc:724
 hloop.cc:725
 hloop.cc:726
 hloop.cc:727
 hloop.cc:728
 hloop.cc:729
 hloop.cc:730
 hloop.cc:731
 hloop.cc:732
 hloop.cc:733
 hloop.cc:734
 hloop.cc:735
 hloop.cc:736
 hloop.cc:737
 hloop.cc:738
 hloop.cc:739
 hloop.cc:740
 hloop.cc:741
 hloop.cc:742
 hloop.cc:743
 hloop.cc:744
 hloop.cc:745
 hloop.cc:746
 hloop.cc:747
 hloop.cc:748
 hloop.cc:749
 hloop.cc:750
 hloop.cc:751
 hloop.cc:752
 hloop.cc:753
 hloop.cc:754
 hloop.cc:755
 hloop.cc:756
 hloop.cc:757
 hloop.cc:758
 hloop.cc:759
 hloop.cc:760
 hloop.cc:761
 hloop.cc:762
 hloop.cc:763
 hloop.cc:764
 hloop.cc:765
 hloop.cc:766
 hloop.cc:767
 hloop.cc:768
 hloop.cc:769
 hloop.cc:770
 hloop.cc:771
 hloop.cc:772
 hloop.cc:773
 hloop.cc:774
 hloop.cc:775
 hloop.cc:776
 hloop.cc:777
 hloop.cc:778
 hloop.cc:779
 hloop.cc:780
 hloop.cc:781
 hloop.cc:782
 hloop.cc:783
 hloop.cc:784
 hloop.cc:785
 hloop.cc:786
 hloop.cc:787
 hloop.cc:788
 hloop.cc:789
 hloop.cc:790
 hloop.cc:791
 hloop.cc:792
 hloop.cc:793
 hloop.cc:794
 hloop.cc:795
 hloop.cc:796
 hloop.cc:797
 hloop.cc:798
 hloop.cc:799
 hloop.cc:800
 hloop.cc:801
 hloop.cc:802
 hloop.cc:803
 hloop.cc:804
 hloop.cc:805
 hloop.cc:806
 hloop.cc:807
 hloop.cc:808
 hloop.cc:809
 hloop.cc:810
 hloop.cc:811
 hloop.cc:812
 hloop.cc:813
 hloop.cc:814
 hloop.cc:815
 hloop.cc:816
 hloop.cc:817
 hloop.cc:818
 hloop.cc:819
 hloop.cc:820
 hloop.cc:821
 hloop.cc:822
 hloop.cc:823
 hloop.cc:824
 hloop.cc:825
 hloop.cc:826
 hloop.cc:827
 hloop.cc:828
 hloop.cc:829
 hloop.cc:830
 hloop.cc:831
 hloop.cc:832
 hloop.cc:833
 hloop.cc:834
 hloop.cc:835
 hloop.cc:836
 hloop.cc:837
 hloop.cc:838
 hloop.cc:839
 hloop.cc:840
 hloop.cc:841
 hloop.cc:842
 hloop.cc:843
 hloop.cc:844
 hloop.cc:845
 hloop.cc:846
 hloop.cc:847
 hloop.cc:848
 hloop.cc:849
 hloop.cc:850
 hloop.cc:851
 hloop.cc:852
 hloop.cc:853
 hloop.cc:854
 hloop.cc:855
 hloop.cc:856
 hloop.cc:857
 hloop.cc:858
 hloop.cc:859
 hloop.cc:860
 hloop.cc:861
 hloop.cc:862
 hloop.cc:863
 hloop.cc:864
 hloop.cc:865
 hloop.cc:866
 hloop.cc:867
 hloop.cc:868
 hloop.cc:869
 hloop.cc:870
 hloop.cc:871
 hloop.cc:872
 hloop.cc:873
 hloop.cc:874
 hloop.cc:875
 hloop.cc:876
 hloop.cc:877
 hloop.cc:878
 hloop.cc:879
 hloop.cc:880
 hloop.cc:881
 hloop.cc:882
 hloop.cc:883
 hloop.cc:884
 hloop.cc:885
 hloop.cc:886
 hloop.cc:887
 hloop.cc:888
 hloop.cc:889
 hloop.cc:890
 hloop.cc:891
 hloop.cc:892
 hloop.cc:893
 hloop.cc:894
 hloop.cc:895
 hloop.cc:896
 hloop.cc:897
 hloop.cc:898
 hloop.cc:899
 hloop.cc:900
 hloop.cc:901
 hloop.cc:902
 hloop.cc:903
 hloop.cc:904
 hloop.cc:905
 hloop.cc:906
 hloop.cc:907
 hloop.cc:908
 hloop.cc:909
 hloop.cc:910
 hloop.cc:911
 hloop.cc:912
 hloop.cc:913
 hloop.cc:914
 hloop.cc:915
 hloop.cc:916
 hloop.cc:917
 hloop.cc:918
 hloop.cc:919
 hloop.cc:920
 hloop.cc:921
 hloop.cc:922
 hloop.cc:923
 hloop.cc:924
 hloop.cc:925
 hloop.cc:926
 hloop.cc:927
 hloop.cc:928
 hloop.cc:929
 hloop.cc:930
 hloop.cc:931
 hloop.cc:932
 hloop.cc:933
 hloop.cc:934
 hloop.cc:935
 hloop.cc:936
 hloop.cc:937
 hloop.cc:938
 hloop.cc:939
 hloop.cc:940
 hloop.cc:941
 hloop.cc:942
 hloop.cc:943
 hloop.cc:944
 hloop.cc:945
 hloop.cc:946
 hloop.cc:947
 hloop.cc:948
 hloop.cc:949
 hloop.cc:950
 hloop.cc:951
 hloop.cc:952
 hloop.cc:953
 hloop.cc:954
 hloop.cc:955
 hloop.cc:956
 hloop.cc:957
 hloop.cc:958
 hloop.cc:959
 hloop.cc:960
 hloop.cc:961
 hloop.cc:962
 hloop.cc:963
 hloop.cc:964
 hloop.cc:965
 hloop.cc:966
 hloop.cc:967
 hloop.cc:968
 hloop.cc:969
 hloop.cc:970
 hloop.cc:971
 hloop.cc:972
 hloop.cc:973
 hloop.cc:974
 hloop.cc:975
 hloop.cc:976
 hloop.cc:977
 hloop.cc:978
 hloop.cc:979
 hloop.cc:980
 hloop.cc:981
 hloop.cc:982
 hloop.cc:983
 hloop.cc:984
 hloop.cc:985
 hloop.cc:986
 hloop.cc:987
 hloop.cc:988
 hloop.cc:989
 hloop.cc:990
 hloop.cc:991
 hloop.cc:992
 hloop.cc:993
 hloop.cc:994
 hloop.cc:995
 hloop.cc:996
 hloop.cc:997
 hloop.cc:998
 hloop.cc:999
 hloop.cc:1000
 hloop.cc:1001
 hloop.cc:1002
 hloop.cc:1003
 hloop.cc:1004
 hloop.cc:1005
 hloop.cc:1006
 hloop.cc:1007
 hloop.cc:1008
 hloop.cc:1009
 hloop.cc:1010
 hloop.cc:1011
 hloop.cc:1012
 hloop.cc:1013
 hloop.cc:1014
 hloop.cc:1015
 hloop.cc:1016
 hloop.cc:1017
 hloop.cc:1018
 hloop.cc:1019
 hloop.cc:1020
 hloop.cc:1021
 hloop.cc:1022
 hloop.cc:1023
 hloop.cc:1024
 hloop.cc:1025
 hloop.cc:1026
 hloop.cc:1027
 hloop.cc:1028
 hloop.cc:1029
 hloop.cc:1030
 hloop.cc:1031
 hloop.cc:1032
 hloop.cc:1033
 hloop.cc:1034
 hloop.cc:1035
 hloop.cc:1036
 hloop.cc:1037
 hloop.cc:1038
 hloop.cc:1039
 hloop.cc:1040
 hloop.cc:1041
 hloop.cc:1042
 hloop.cc:1043