00001
00002
00003
00004
00005
00006
00007 #include "Stiostream.h"
00008 #include "St_sfs_Maker.h"
00009 #include "StHit.h"
00010 #include "StEventTypes.h"
00011 #include "tables/St_HitError_Table.h"
00012 #include "StSvtHitCollection.h"
00013 #include "StEvent.h"
00014 #include "TMath.h"
00015 #include "TRandom.h"
00016 #include "StSvtDbMaker/StSvtDbMaker.h"
00017 #include "tables/St_g2t_svt_hit_Table.h"
00018 #include "TGeoMatrix.h"
00019 ClassImp(St_sfs_Maker);
00020
00021 Int_t St_sfs_Maker::InitRun(Int_t RunNo)
00022 {
00023
00024 St_HitError *tableSet = (St_HitError *) GetDataBase("Calibrations/tracker/svtHitError");
00025 HitError_st* hitError = tableSet->GetTable();
00026 mResXSvt = TMath::Sqrt(hitError->coeff[0]);
00027 mResZSvt = TMath::Sqrt(hitError->coeff[3]);
00028 LOG_DEBUG << "Smearing SVT hits by " << mResXSvt << " " << mResZSvt << " cm in the SVT " << endm;
00029 return kStOk;
00030 }
00031
00032 Int_t St_sfs_Maker::Make()
00033 {
00034 if (! gRandom) gRandom = new TRandom();
00035
00036 StEvent *rEvent = (StEvent*) GetInputDS("StEvent");
00037 if (! rEvent) { LOG_WARN << "No StEvent on input, bye bye" << endm; return kStWarn; }
00038 StSvtHitCollection *rCol = rEvent->svtHitCollection();
00039 if (!rCol) {
00040 rCol = new StSvtHitCollection;
00041 rEvent->setSvtHitCollection(rCol);
00042 }
00043 TDataSetIter geant(GetInputDS("geant"));
00044 St_g2t_svt_hit *g2t_svt_hit = (St_g2t_svt_hit *) geant("g2t_svt_hit");
00045 if (! g2t_svt_hit) {
00046 LOG_WARN << "No g2t_svt_hit on input, bye bye" << endm; return kStWarn;
00047 return kStWarn;
00048 }
00049 g2t_svt_hit_st *g2t = g2t_svt_hit->GetTable();
00050 Int_t Nhits = g2t_svt_hit->GetNRows();
00051 if (Nhits <= 0) return kStWarn;
00052 THashList *rotMHash = gStSvtDbMaker->GetRotations();
00053 if (! rotMHash) {LOG_WARN << " No list of Rotation from StSvtDbMaker " << endm; return kStWarn;}
00054
00055 for (Int_t i = 0; i < Nhits; i++) {
00056 TGeoHMatrix *comb = (TGeoHMatrix *) rotMHash->FindObject(Form("R%04i",g2t[i].volume_id%10000));
00057 if (! comb) { LOG_WARN << "No rotation matrix for wafer =" << g2t[i].volume_id << endm; continue;}
00058 Double_t xyzG[3] = {g2t[i].x[0],g2t[i].x[1],g2t[i].x[2]};
00059 Double_t xyzL[3];
00060 comb->MasterToLocal(xyzG,xyzL);
00061 if (Debug() && TMath::Abs(xyzL[2]) > 0.0100) {
00062 cout << "Id: " << g2t[i].volume_id
00063 << "\txyzL :" << xyzL[0] << "\t" << xyzL[1] << "\t" << xyzL[2] << endl;
00064 comb->Print();
00065 }
00066 xyzL[0] += gRandom->Gaus(0, mResXSvt);
00067 xyzL[1] += gRandom->Gaus(0, mResZSvt);
00068 comb->LocalToMaster(xyzL,xyzG);
00069
00070 Int_t ladderID = g2t[i].volume_id%100;
00071 Int_t waferID = (g2t[i].volume_id/100)%10;
00072 Int_t layerID = (g2t[i].volume_id%10000)/1000;
00073 Int_t barrelID = (layerID-1)/2 + 1;
00074 Int_t hybridID = 1;
00075 if (xyzL[0] >= 0) hybridID = 2;
00076 static const Int_t NumberOfHybrids = 2;
00077 static const Int_t NumberOfWafers[3] = {4, 6, 7};
00078 static const Int_t NumberOfLadders[3] = {8,12,16};
00079 Int_t index = ((ladderID-1)*NumberOfWafers[barrelID-1] + (waferID-1))*NumberOfHybrids + (hybridID-1);
00080 switch (barrelID) {
00081 case 3: index += NumberOfLadders[barrelID-3]*NumberOfWafers[barrelID-3]*NumberOfHybrids;
00082 case 2: index += NumberOfLadders[barrelID-2]*NumberOfWafers[barrelID-2]*NumberOfHybrids;
00083 case 1: break;
00084 default:
00085 LOG_ERROR << "There is NO barrel number " << barrelID << " !!!" << endm;
00086
00087 index = -1;
00088 break;
00089 }
00090 assert(index != -1);
00091 Int_t hw = 2;
00092 hw += (1L<<4)*(index);
00093 StSvtHit *svtHit = new StSvtHit();
00094 svtHit->setHardwarePosition(hw);
00095 svtHit->setLocalPosition(xyzL[0],xyzL[1]);
00096 svtHit->setCharge(g2t[i].de);
00097 svtHit->setIdTruth(g2t[i].track_p,100);
00098 svtHit->setFlag(0);
00099 svtHit->setFitFlag(0);
00100 StThreeVectorF pos(xyzG[0],xyzG[1],xyzG[2]);
00101 svtHit->setPosition(pos);
00102 rCol->addHit(svtHit);
00103 }
00104 return kStOK;
00105 }
