ROOT logo
//////////////////////////////////////////////////////////////////////////////
//
// $Id: $
//
//*-- Author  : S. Lebedev
//
//_HADES_CLASS_DESCRIPTION
//////////////////////////////////////////////////////////////////////////////
//
//  HRich700RecoQa
//
//
//////////////////////////////////////////////////////////////////////////////


#include "hrich700recoqa.h"
#include "hades.h"
#include "hcategory.h"
#include "hevent.h"
#include "hlinearcatiter.h"
#include "hmatrixcatiter.h"
#include "hparset.h"
#include "hspectrometer.h"
#include "richdef.h"
#include "hrichcalsim.h"
#include "hgeantkine.h"
#include "hgeantrich.h"
#include "hrich700histmanager.h"
#include "hrichhitsim.h"
#include "hrich700ringfittercop.h"
#include "hrich700drawhist.h"
#include "hrich700utils.h"



#include "TCanvas.h"

#include <iostream>
#include <sstream>
#include <map>
#include <vector>

using namespace std;


ClassImp(HRich700RecoQa)

    HRich700RecoQa::HRich700RecoQa():
fEventNum(0),
    fMinNofRichCalsAcc(5),
    fTrueQuota(0.7)
{

}

HRich700RecoQa::~HRich700RecoQa()
{

}

Bool_t HRich700RecoQa::init()
{

    fCatKine = gHades->getCurrentEvent()->getCategory(catGeantKine);
    if (NULL == fCatKine) {
	Error("init", "Initializatin of kine category failed, returning...");
	return kFALSE;
    }

    fCatRichHit = gHades->getCurrentEvent()->getCategory(catRichHit);
    if (NULL == fCatRichHit) {
	Error("init", "Initialization of geant category for catRichHit failed, returning...");
	return kFALSE;
    }

    fCatRichCal = gHades->getCurrentEvent()->getCategory(catRichCal);
    if (NULL == fCatRichCal) {
	Error("init", "Initialization of geant category for catRichCal failed, returning...");
	return kFALSE;
    }

    initHist();

    return kTRUE;

}


Bool_t HRich700RecoQa::reinit()
{
    return kTRUE;
}

void HRich700RecoQa::initHist()
{
    fHM = new HRich700HistManager();
    fHM->Create1<TH1D>("fhTrueAllRatioAll"            , "fhTrueAllRatioAll;true/all;Yield"                        , 110, -0.05, 1.05);

    fHM->Create1<TH1D>("fhTrueAllRatioElectron"       , "fhTrueAllRatioElectron;true/all;Yield"                   , 110, -0.05, 1.05);
    fHM->Create2<TH2D>("fhTrueAllRatioVsThetaElectron", "fhTrueAllRatioVsThetaElectron;Theta [deg];true/all;Yield", 30 ,  0.  , 90., 55, -0.05, 1.05);
    fHM->Create1<TH1D>("fhNofTrueCalsElectron"        , "fhNofTrueCalsElectron;Nof cals;Yield"                    , 40 , -0.5 , 39.5);
    fHM->Create1<TH1D>("fhNofWrongCalsElectron"       , "fhNofWrongCalsElectron;Nof cals;Yield"                   , 40 , -0.5 , 39.5);
    fHM->Create1<TH1D>("fhNofAllCalsElectron"         , "fhNofAllCalsElectron;Nof cals;Yield"                     , 40 , -0.5 , 39.5);

    fHM->Create1<TH1D>("fhElectronAccTheta"           , "fhElectronAccTheta;Theta [deg];Yield"                    , 90 ,  0.  , 90.);
    fHM->Create1<TH1D>("fhElectronRecTheta"           , "fhElectronRecTheta;Theta [deg];Yield"                    , 90 ,  0.  , 90.);
    fHM->Create1<TH1D>("fhElectronAccPhi"             , "fhElectronAccPhi;Phi [deg];Yield"                        , 90 ,  0.  , 360.);
    fHM->Create1<TH1D>("fhElectronRecPhi"             , "fhElectronRecPhi;Phi [deg];Yield"                        , 90 ,  0.  , 360.);
    fHM->Create1<TH1D>("fhElectronAccMom"             , "fhElectronAccMom;P [MeV/c];Yield"                        , 100,  0.  , 1500.);
    fHM->Create1<TH1D>("fhElectronRecMom"             , "fhElectronRecMom;P [MeV/c];Yield"                        , 100,  0.  , 1500.);
    fHM->Create1<TH1D>("fhElectronAccNofCals"         , "fhElectronAccNofCals;Number of cals;Yield"               , 30 , -0.5 , 29.5);
    fHM->Create1<TH1D>("fhElectronRecNofCals"         , "fhElectronRecNofCals;Number of cals;Yield"               , 30 , -0.5 , 29.5);

    fHM->Create1<TH1D>("fhGhostsClones"               , "fhGhostsClones;Ghosts and clones;Number per event"       , 2  , -0.5 , 1.5);

    //pair efficiency
    fHM->Create1<TH1D>("fhPairAccTheta"               , "fhPairAccTheta;Theta [deg];Yield"                        , 90 , 0.   , 90.);
    fHM->Create1<TH1D>("fhPairRecTheta"               , "fhPairRecTheta;Theta [deg];Yield"                        , 90 , 0.   , 90.);
}

void HRich700RecoQa::initRichHitMap()
{
    fHitMap.clear();
    Int_t nofRichCals = fCatRichCal->getEntries();
    for (Int_t iC = 0; iC < nofRichCals; iC++) {
	HRichCalSim* calSim = static_cast<HRichCalSim*>(fCatRichCal->getObject(iC));
	if (NULL == calSim) continue;
	Int_t nofTrackIds = calSim->getNofTrackIds();
	for (Int_t iT = 0; iT < nofTrackIds; iT++) {
	    Int_t trackId = calSim->getTrackId(iT);
	    if (trackId < 0) continue;
	    fHitMap[trackId]++;
	}
    }
}


Int_t HRich700RecoQa::execute()
{
    fEventNum++;
    initRichHitMap();
    processEvent();
    fillAccRecHist();

    return 0;
}

void HRich700RecoQa::processEvent()
{
    Int_t nofRichRings = fCatRichHit->getEntries();
    for (Int_t i = 0; i < nofRichRings; i++) {
	HRichHitSim* ring = static_cast<HRichHitSim*>(fCatRichHit->getObject(i));
	if (ring == NULL) continue;
	Int_t trueCals = ring->weigTrack1;
	Int_t trackId = ring->track1;
	Int_t nofCals = ring->fRich700NofRichCals;
	Double_t trueAllRatio = (Double_t)trueCals / (Double_t) nofCals;

	if (trackId <= 0) continue;
	fHM->H1("fhTrueAllRatioAll")->Fill(trueAllRatio);
	HGeantKine* kine = (HGeantKine*)fCatKine-> getObject(trackId - 1);
	if (isPrimaryElectron(kine)) {
	    fHM->H1("fhTrueAllRatioElectron")->Fill(trueAllRatio);
	    fHM->H2("fhTrueAllRatioVsThetaElectron")->Fill(kine->getThetaDeg(), trueAllRatio);
	    fHM->H1("fhNofTrueCalsElectron")->Fill(trueCals);
	    fHM->H1("fhNofWrongCalsElectron")->Fill(nofCals - trueCals);
	    fHM->H1("fhNofAllCalsElectron")->Fill(nofCals);
	}
    }
}

void HRich700RecoQa::fillAccRecHist()
{
    Int_t nofAccEl = 0;
    Double_t thetaAcc[2];
    Int_t nofKine = fCatKine->getEntries();
    for (Int_t i = 0; i < nofKine; i++) {
	HGeantKine* kine = (HGeantKine*)fCatKine-> getObject(i);
	Double_t theta = kine->getThetaDeg();
	Double_t phi = kine->getPhiDeg();
	Double_t mom = kine->getTotalMomentum();
	Int_t nofCals = fHitMap[kine->getTrack()];
	if (isPrimaryElectron(kine) && isRichAcc(kine->getTrack())) {
	    fHM->H1("fhElectronAccTheta")->Fill(theta);
	    fHM->H1("fhElectronAccPhi")->Fill(phi);
	    fHM->H1("fhElectronAccMom")->Fill(mom);
	    fHM->H1("fhElectronAccNofCals")->Fill(nofCals);
	    if (nofAccEl <= 1) thetaAcc[nofAccEl] = theta;
	    nofAccEl++;
	}
    }
    if (nofAccEl == 2){
	fHM->H1("fhPairAccTheta")->Fill(0.5 * (thetaAcc[0] + thetaAcc[1]));
    }

    Int_t nofRecEl = 0;
    Double_t thetaRec[2];
    Int_t nofRichRings = fCatRichHit->getEntries();
    map<Int_t, Int_t> clonesMap;
    for (Int_t i = 0; i < nofRichRings; i++) {
	HRichHitSim* ring = static_cast<HRichHitSim*>(fCatRichHit->getObject(i));
	if (ring == NULL) continue;
	Int_t trackId = ring->track1;
	Double_t trueAllRatio = (Double_t)ring->weigTrack1 / (Double_t) ring->fRich700NofRichCals;
	Bool_t isTrue = (trueAllRatio >= fTrueQuota);
	if (trackId <= 0){
	    fHM->H1("fhGhostsClones")->Fill(0.);
	    continue;
	}
	HGeantKine* kine = (HGeantKine*)fCatKine-> getObject(trackId - 1);
	if (kine == NULL) continue;
	Double_t theta = kine->getThetaDeg();
	Double_t phi = kine->getPhiDeg();
	Double_t mom = kine->getTotalMomentum();
	Int_t nofCals = fHitMap[trackId];
	Bool_t isClone = (clonesMap[trackId] > 0);
	if (isPrimaryElectron(kine) && isRichAcc(trackId) && isTrue && !isClone) {
	    fHM->H1("fhElectronRecTheta")->Fill(theta);
	    fHM->H1("fhElectronRecPhi")->Fill(phi);
	    fHM->H1("fhElectronRecMom")->Fill(mom);
	    fHM->H1("fhElectronRecNofCals")->Fill(nofCals);
	    if (nofRecEl <= 1) thetaRec[nofRecEl] = theta;
	    nofRecEl++;
	    clonesMap[trackId]++;
	}

	if (!isTrue) {
	    fHM->H1("fhGhostsClones")->Fill(0.);
	}
	if (isClone) {
	    fHM->H1("fhGhostsClones")->Fill(1.);
	}
    }
    if (nofRecEl == 2){
	fHM->H1("fhPairRecTheta")->Fill(0.5*(thetaRec[0] + thetaRec[1]));
    }
}

Bool_t HRich700RecoQa::isRichAcc(
				 Int_t trackId)
{
    if (fHitMap[trackId] >= fMinNofRichCalsAcc) return kTRUE;
    return kFALSE;
}

Bool_t HRich700RecoQa::isPrimaryElectron(
					 HGeantKine* kine)
{
    if (kine == NULL) return kFALSE;
    return ((kine->getID() == 2 || kine->getID() == 3) && kine->isPrimary());
}

void HRich700RecoQa::drawHist()
{
    HRichDrawHist::SetDefaultDrawStyle();
    {
	TCanvas* c = fHM->CreateCanvas("hrich_fhTrueAllRatioAll", "hrich_fhTrueAllRatioAll", 800, 800);
	c->cd();
	fHM->NormalizeToIntegral("fhTrueAllRatioAll");
	HRichDrawHist::DrawH1(fHM->H1("fhTrueAllRatioAll"));
    }

    {
	TCanvas* c = fHM->CreateCanvas("hrich_fhTrueAllRatioElectron", "hrich_fhTrueAllRatioElectron", 800, 800);
	c->cd();
	fHM->NormalizeToIntegral("fhTrueAllRatioElectron");
	HRichDrawHist::DrawH1(fHM->H1("fhTrueAllRatioElectron"));
    }

    {
	TCanvas* c = fHM->CreateCanvas("hrich_fhTrueAllRatioVsThetaElectron", "hrich_fhTrueAllRatioVsThetaElectron", 1000, 800);
	c->cd();
	fHM->NormalizeToIntegral("fhTrueAllRatioVsThetaElectron");
	HRichDrawHist::DrawH2(fHM->H2("fhTrueAllRatioVsThetaElectron"));
    }

    {
	TCanvas* c = fHM->CreateCanvas("hrich_fhNofCalsElectron", "hrich_fhNofCalsElectron", 800, 800);
	c->cd();
	fHM->NormalizeToIntegral("fhNofTrueCalsElectron");
	fHM->NormalizeToIntegral("fhNofWrongCalsElectron");
	fHM->NormalizeToIntegral("fhNofAllCalsElectron");
	vector<TH1*> hist;
	hist.push_back(fHM->H1("fhNofTrueCalsElectron"));
	hist.push_back(fHM->H1("fhNofWrongCalsElectron"));
	hist.push_back(fHM->H1("fhNofAllCalsElectron"));
	vector<string> labels;
	stringstream ss;
	ss.precision(2);
	ss << "True (" << fixed << fHM->H1("fhNofTrueCalsElectron")->GetMean() << ")";
	labels.push_back(ss.str());
	ss.str("");
	ss << "Wrong (" << fixed << fHM->H1("fhNofWrongCalsElectron")->GetMean() << ")";
	labels.push_back(ss.str());
	ss.str("");
	ss << "All (" << fixed << fHM->H1("fhNofAllCalsElectron")->GetMean() << ")";
	labels.push_back(ss.str());
	HRichDrawHist::DrawH1(hist, labels);
	gPad->SetLogy(kTRUE);
    }

    drawAccRecEff("hrich_fhElectronAccRecTheta", "hrich_fhElectronEffTheta", "fhElectronAccTheta", "fhElectronRecTheta");
    drawAccRecEff("hrich_fhElectronAccRecPhi", "hrich_fhElectronEffPhi", "fhElectronAccPhi", "fhElectronRecPhi");
    drawAccRecEff("hrich_fhElectronAccRecMom", "hrich_fhElectronEffMom", "fhElectronAccMom", "fhElectronRecMom");
    drawAccRecEff("hrich_fhElectronAccRecNofCals", "hrich_fhElectronEffNofCals", "fhElectronAccNofCals", "fhElectronRecNofCals");

    drawAccRecEff("hrich_fhPairAccRecTheta", "hrich_fhPairEffTheta", "fhPairAccTheta", "fhPairRecTheta");

    {
	TCanvas* c = fHM->CreateCanvas("hrich_fhGhostsClones", "hrich_fhGhostsClones", 800, 800);
	c->cd();
	fHM->Scale("fhGhostsClones", 1./ (Double_t)fEventNum);
	HRichDrawHist::DrawH1(fHM->H1("fhGhostsClones"));
    }

}

void HRich700RecoQa::drawAccRecEff(
				   const string& canvasNameAccRec,
				   const string& canvasNameEff,
				   const string& histNameAcc,
				   const string& histNameRec)
{
    {
	TCanvas* c = fHM->CreateCanvas(canvasNameAccRec.c_str(), canvasNameAccRec.c_str(), 800, 800);
	c->cd();
	vector<TH1*> hist;
	hist.push_back(fHM->H1(histNameAcc));
	hist.push_back(fHM->H1(histNameRec));
	vector<string> labels;
	stringstream ss;
	ss.precision(3);
	ss << "Acc (" << fixed << fHM->H1(histNameAcc)->GetEntries() / (Double_t)fEventNum << ")";
	labels.push_back(ss.str());
	ss.str("");
	ss << "Rec (" << fixed << fHM->H1(histNameRec)->GetEntries() / (Double_t)fEventNum << ")";
	labels.push_back(ss.str());
	HRichDrawHist::DrawH1(hist, labels);
    }

    {
	TH1D* pxEff = RichUtils::DivideH1((TH1D*)fHM->H1(histNameRec)->Clone(), (TH1D*)fHM->H1(histNameAcc)->Clone(), "", 100., "Efficiency [%]");
	pxEff->SetMinimum(0.);
	TCanvas* c = fHM->CreateCanvas(canvasNameEff.c_str(), canvasNameEff.c_str(), 800, 800);
	c->cd();
	vector<TH1*> hist;
	hist.push_back(pxEff);
	vector<string> labels;
	stringstream ss;
	ss.precision(1);
	ss << "Electrons (" << fixed << 100.*fHM->H1(histNameRec)->GetEntries() / fHM->H1(histNameAcc)->GetEntries() << "%)";
	labels.push_back(ss.str());
	HRichDrawHist::DrawH1(hist, labels);
    }
}

Bool_t HRich700RecoQa::finalize()
{
    drawHist();
    fHM->SaveCanvasToImage(string(fOutputDir + "/"));
    return kTRUE;
}
 hrich700recoqa.cc:1
 hrich700recoqa.cc:2
 hrich700recoqa.cc:3
 hrich700recoqa.cc:4
 hrich700recoqa.cc:5
 hrich700recoqa.cc:6
 hrich700recoqa.cc:7
 hrich700recoqa.cc:8
 hrich700recoqa.cc:9
 hrich700recoqa.cc:10
 hrich700recoqa.cc:11
 hrich700recoqa.cc:12
 hrich700recoqa.cc:13
 hrich700recoqa.cc:14
 hrich700recoqa.cc:15
 hrich700recoqa.cc:16
 hrich700recoqa.cc:17
 hrich700recoqa.cc:18
 hrich700recoqa.cc:19
 hrich700recoqa.cc:20
 hrich700recoqa.cc:21
 hrich700recoqa.cc:22
 hrich700recoqa.cc:23
 hrich700recoqa.cc:24
 hrich700recoqa.cc:25
 hrich700recoqa.cc:26
 hrich700recoqa.cc:27
 hrich700recoqa.cc:28
 hrich700recoqa.cc:29
 hrich700recoqa.cc:30
 hrich700recoqa.cc:31
 hrich700recoqa.cc:32
 hrich700recoqa.cc:33
 hrich700recoqa.cc:34
 hrich700recoqa.cc:35
 hrich700recoqa.cc:36
 hrich700recoqa.cc:37
 hrich700recoqa.cc:38
 hrich700recoqa.cc:39
 hrich700recoqa.cc:40
 hrich700recoqa.cc:41
 hrich700recoqa.cc:42
 hrich700recoqa.cc:43
 hrich700recoqa.cc:44
 hrich700recoqa.cc:45
 hrich700recoqa.cc:46
 hrich700recoqa.cc:47
 hrich700recoqa.cc:48
 hrich700recoqa.cc:49
 hrich700recoqa.cc:50
 hrich700recoqa.cc:51
 hrich700recoqa.cc:52
 hrich700recoqa.cc:53
 hrich700recoqa.cc:54
 hrich700recoqa.cc:55
 hrich700recoqa.cc:56
 hrich700recoqa.cc:57
 hrich700recoqa.cc:58
 hrich700recoqa.cc:59
 hrich700recoqa.cc:60
 hrich700recoqa.cc:61
 hrich700recoqa.cc:62
 hrich700recoqa.cc:63
 hrich700recoqa.cc:64
 hrich700recoqa.cc:65
 hrich700recoqa.cc:66
 hrich700recoqa.cc:67
 hrich700recoqa.cc:68
 hrich700recoqa.cc:69
 hrich700recoqa.cc:70
 hrich700recoqa.cc:71
 hrich700recoqa.cc:72
 hrich700recoqa.cc:73
 hrich700recoqa.cc:74
 hrich700recoqa.cc:75
 hrich700recoqa.cc:76
 hrich700recoqa.cc:77
 hrich700recoqa.cc:78
 hrich700recoqa.cc:79
 hrich700recoqa.cc:80
 hrich700recoqa.cc:81
 hrich700recoqa.cc:82
 hrich700recoqa.cc:83
 hrich700recoqa.cc:84
 hrich700recoqa.cc:85
 hrich700recoqa.cc:86
 hrich700recoqa.cc:87
 hrich700recoqa.cc:88
 hrich700recoqa.cc:89
 hrich700recoqa.cc:90
 hrich700recoqa.cc:91
 hrich700recoqa.cc:92
 hrich700recoqa.cc:93
 hrich700recoqa.cc:94
 hrich700recoqa.cc:95
 hrich700recoqa.cc:96
 hrich700recoqa.cc:97
 hrich700recoqa.cc:98
 hrich700recoqa.cc:99
 hrich700recoqa.cc:100
 hrich700recoqa.cc:101
 hrich700recoqa.cc:102
 hrich700recoqa.cc:103
 hrich700recoqa.cc:104
 hrich700recoqa.cc:105
 hrich700recoqa.cc:106
 hrich700recoqa.cc:107
 hrich700recoqa.cc:108
 hrich700recoqa.cc:109
 hrich700recoqa.cc:110
 hrich700recoqa.cc:111
 hrich700recoqa.cc:112
 hrich700recoqa.cc:113
 hrich700recoqa.cc:114
 hrich700recoqa.cc:115
 hrich700recoqa.cc:116
 hrich700recoqa.cc:117
 hrich700recoqa.cc:118
 hrich700recoqa.cc:119
 hrich700recoqa.cc:120
 hrich700recoqa.cc:121
 hrich700recoqa.cc:122
 hrich700recoqa.cc:123
 hrich700recoqa.cc:124
 hrich700recoqa.cc:125
 hrich700recoqa.cc:126
 hrich700recoqa.cc:127
 hrich700recoqa.cc:128
 hrich700recoqa.cc:129
 hrich700recoqa.cc:130
 hrich700recoqa.cc:131
 hrich700recoqa.cc:132
 hrich700recoqa.cc:133
 hrich700recoqa.cc:134
 hrich700recoqa.cc:135
 hrich700recoqa.cc:136
 hrich700recoqa.cc:137
 hrich700recoqa.cc:138
 hrich700recoqa.cc:139
 hrich700recoqa.cc:140
 hrich700recoqa.cc:141
 hrich700recoqa.cc:142
 hrich700recoqa.cc:143
 hrich700recoqa.cc:144
 hrich700recoqa.cc:145
 hrich700recoqa.cc:146
 hrich700recoqa.cc:147
 hrich700recoqa.cc:148
 hrich700recoqa.cc:149
 hrich700recoqa.cc:150
 hrich700recoqa.cc:151
 hrich700recoqa.cc:152
 hrich700recoqa.cc:153
 hrich700recoqa.cc:154
 hrich700recoqa.cc:155
 hrich700recoqa.cc:156
 hrich700recoqa.cc:157
 hrich700recoqa.cc:158
 hrich700recoqa.cc:159
 hrich700recoqa.cc:160
 hrich700recoqa.cc:161
 hrich700recoqa.cc:162
 hrich700recoqa.cc:163
 hrich700recoqa.cc:164
 hrich700recoqa.cc:165
 hrich700recoqa.cc:166
 hrich700recoqa.cc:167
 hrich700recoqa.cc:168
 hrich700recoqa.cc:169
 hrich700recoqa.cc:170
 hrich700recoqa.cc:171
 hrich700recoqa.cc:172
 hrich700recoqa.cc:173
 hrich700recoqa.cc:174
 hrich700recoqa.cc:175
 hrich700recoqa.cc:176
 hrich700recoqa.cc:177
 hrich700recoqa.cc:178
 hrich700recoqa.cc:179
 hrich700recoqa.cc:180
 hrich700recoqa.cc:181
 hrich700recoqa.cc:182
 hrich700recoqa.cc:183
 hrich700recoqa.cc:184
 hrich700recoqa.cc:185
 hrich700recoqa.cc:186
 hrich700recoqa.cc:187
 hrich700recoqa.cc:188
 hrich700recoqa.cc:189
 hrich700recoqa.cc:190
 hrich700recoqa.cc:191
 hrich700recoqa.cc:192
 hrich700recoqa.cc:193
 hrich700recoqa.cc:194
 hrich700recoqa.cc:195
 hrich700recoqa.cc:196
 hrich700recoqa.cc:197
 hrich700recoqa.cc:198
 hrich700recoqa.cc:199
 hrich700recoqa.cc:200
 hrich700recoqa.cc:201
 hrich700recoqa.cc:202
 hrich700recoqa.cc:203
 hrich700recoqa.cc:204
 hrich700recoqa.cc:205
 hrich700recoqa.cc:206
 hrich700recoqa.cc:207
 hrich700recoqa.cc:208
 hrich700recoqa.cc:209
 hrich700recoqa.cc:210
 hrich700recoqa.cc:211
 hrich700recoqa.cc:212
 hrich700recoqa.cc:213
 hrich700recoqa.cc:214
 hrich700recoqa.cc:215
 hrich700recoqa.cc:216
 hrich700recoqa.cc:217
 hrich700recoqa.cc:218
 hrich700recoqa.cc:219
 hrich700recoqa.cc:220
 hrich700recoqa.cc:221
 hrich700recoqa.cc:222
 hrich700recoqa.cc:223
 hrich700recoqa.cc:224
 hrich700recoqa.cc:225
 hrich700recoqa.cc:226
 hrich700recoqa.cc:227
 hrich700recoqa.cc:228
 hrich700recoqa.cc:229
 hrich700recoqa.cc:230
 hrich700recoqa.cc:231
 hrich700recoqa.cc:232
 hrich700recoqa.cc:233
 hrich700recoqa.cc:234
 hrich700recoqa.cc:235
 hrich700recoqa.cc:236
 hrich700recoqa.cc:237
 hrich700recoqa.cc:238
 hrich700recoqa.cc:239
 hrich700recoqa.cc:240
 hrich700recoqa.cc:241
 hrich700recoqa.cc:242
 hrich700recoqa.cc:243
 hrich700recoqa.cc:244
 hrich700recoqa.cc:245
 hrich700recoqa.cc:246
 hrich700recoqa.cc:247
 hrich700recoqa.cc:248
 hrich700recoqa.cc:249
 hrich700recoqa.cc:250
 hrich700recoqa.cc:251
 hrich700recoqa.cc:252
 hrich700recoqa.cc:253
 hrich700recoqa.cc:254
 hrich700recoqa.cc:255
 hrich700recoqa.cc:256
 hrich700recoqa.cc:257
 hrich700recoqa.cc:258
 hrich700recoqa.cc:259
 hrich700recoqa.cc:260
 hrich700recoqa.cc:261
 hrich700recoqa.cc:262
 hrich700recoqa.cc:263
 hrich700recoqa.cc:264
 hrich700recoqa.cc:265
 hrich700recoqa.cc:266
 hrich700recoqa.cc:267
 hrich700recoqa.cc:268
 hrich700recoqa.cc:269
 hrich700recoqa.cc:270
 hrich700recoqa.cc:271
 hrich700recoqa.cc:272
 hrich700recoqa.cc:273
 hrich700recoqa.cc:274
 hrich700recoqa.cc:275
 hrich700recoqa.cc:276
 hrich700recoqa.cc:277
 hrich700recoqa.cc:278
 hrich700recoqa.cc:279
 hrich700recoqa.cc:280
 hrich700recoqa.cc:281
 hrich700recoqa.cc:282
 hrich700recoqa.cc:283
 hrich700recoqa.cc:284
 hrich700recoqa.cc:285
 hrich700recoqa.cc:286
 hrich700recoqa.cc:287
 hrich700recoqa.cc:288
 hrich700recoqa.cc:289
 hrich700recoqa.cc:290
 hrich700recoqa.cc:291
 hrich700recoqa.cc:292
 hrich700recoqa.cc:293
 hrich700recoqa.cc:294
 hrich700recoqa.cc:295
 hrich700recoqa.cc:296
 hrich700recoqa.cc:297
 hrich700recoqa.cc:298
 hrich700recoqa.cc:299
 hrich700recoqa.cc:300
 hrich700recoqa.cc:301
 hrich700recoqa.cc:302
 hrich700recoqa.cc:303
 hrich700recoqa.cc:304
 hrich700recoqa.cc:305
 hrich700recoqa.cc:306
 hrich700recoqa.cc:307
 hrich700recoqa.cc:308
 hrich700recoqa.cc:309
 hrich700recoqa.cc:310
 hrich700recoqa.cc:311
 hrich700recoqa.cc:312
 hrich700recoqa.cc:313
 hrich700recoqa.cc:314
 hrich700recoqa.cc:315
 hrich700recoqa.cc:316
 hrich700recoqa.cc:317
 hrich700recoqa.cc:318
 hrich700recoqa.cc:319
 hrich700recoqa.cc:320
 hrich700recoqa.cc:321
 hrich700recoqa.cc:322
 hrich700recoqa.cc:323
 hrich700recoqa.cc:324
 hrich700recoqa.cc:325
 hrich700recoqa.cc:326
 hrich700recoqa.cc:327
 hrich700recoqa.cc:328
 hrich700recoqa.cc:329
 hrich700recoqa.cc:330
 hrich700recoqa.cc:331
 hrich700recoqa.cc:332
 hrich700recoqa.cc:333
 hrich700recoqa.cc:334
 hrich700recoqa.cc:335
 hrich700recoqa.cc:336
 hrich700recoqa.cc:337
 hrich700recoqa.cc:338
 hrich700recoqa.cc:339
 hrich700recoqa.cc:340
 hrich700recoqa.cc:341
 hrich700recoqa.cc:342
 hrich700recoqa.cc:343
 hrich700recoqa.cc:344
 hrich700recoqa.cc:345
 hrich700recoqa.cc:346
 hrich700recoqa.cc:347
 hrich700recoqa.cc:348
 hrich700recoqa.cc:349
 hrich700recoqa.cc:350
 hrich700recoqa.cc:351
 hrich700recoqa.cc:352
 hrich700recoqa.cc:353
 hrich700recoqa.cc:354
 hrich700recoqa.cc:355
 hrich700recoqa.cc:356
 hrich700recoqa.cc:357
 hrich700recoqa.cc:358
 hrich700recoqa.cc:359
 hrich700recoqa.cc:360
 hrich700recoqa.cc:361
 hrich700recoqa.cc:362