Newer
Older
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
MyMainFrame::MyMainFrame(const Setup* config_, const TGWindow *p,UInt_t w,UInt_t h): config(config_) {
// use hierarchical cleaning
fMain = new TGMainFrame(p,w,h);
fMain->SetCleanup(kDeepCleanup);
fMain->Connect("CloseWindow()", "MyMainFrame", this, "CloseWindow()");
// create layout
Init();
canvas = fEcanvas->GetCanvas();
canvas->Update();
// popup window
fMain->SetWindowName("Back projection reconstruction");
fMain->MapSubwindows();
fMain->Resize(fMain->GetDefaultSize());
fMain->Resize(w, h);
fMain->MapWindow();
}
MyMainFrame::~MyMainFrame() {
// // Clean up used widgets: frames, buttons, layout hints
// fMain->Cleanup();
//buttons
delete fF3;
delete fL3;
delete fStopB;
delete fClearB;
delete fStartB;
// canvas
delete fF5;
delete fL4;
delete fEcanvas;
delete fMain;
delete histo;
for (int i = 0; i < Nl; i++)
{
delete latitudes[i];
}
for (int i = 0; i < NL; i++)
{
delete longitudes[i];
}
delete sourceMarker;
}
void MyMainFrame::CloseWindow() {
// Got close message for this MainFrame. Terminates the application.
gApplication->Terminate(0);
}
int MyMainFrame::Init(){
// buttons
fF3 = new TGCompositeFrame(fMain, 60, 20, kHorizontalFrame);
fStartB = new TGTextButton(fF3, "Start");
fStartB->Connect("Clicked()", "MyMainFrame", this, "Start()");
fClearB = new TGTextButton(fF3, "Clear");
fClearB->Connect("Clicked()", "MyMainFrame", this, "Clear()");
fStopB = new TGTextButton(fF3, "Stop");
fStopB->Connect("Clicked()", "MyMainFrame", this, "Stop()");
// // counts
// fCountsLabel = new TGLabel(fF3, "Total counts");
// fCounts = new TGNumberEntryField(fF3, 1, 0, TGNumberFormat::kNESInteger,
// TGNumberEntryField::kNEANonNegative);
// fCounts->SetAlignment(kTextLeft);
// fCounts->SetLimits(TGNumberEntryField::kNELLimitMin, 0);
// fCounts->SetEditable(kFALSE);
// // fCounts->Resize(100,20);
fL3 = new TGLayoutHints(kLHintsTop | kLHintsLeft,
5, 5, 5, 5);
fF3->AddFrame(fStartB, fL3);
fF3->AddFrame(fClearB, fL3);
fF3->AddFrame(fStopB, fL3);
// fF3->AddFrame(fCountsLabel, new TGLayoutHints(kLHintsCenterY | kLHintsLeft, 2, 2, 2, 2));
// fF3->AddFrame(fCounts, fL3);
// embedded canvas
fF5 = new TGCompositeFrame(fMain, 60, 60, kHorizontalFrame);
fL4 = new TGLayoutHints(kLHintsTop | kLHintsLeft
| kLHintsExpandX | kLHintsExpandY, 5, 5, 5, 5);
fEcanvas = new TRootEmbeddedCanvas("ec1", fF5, 100, 100);
fF5->AddFrame(fEcanvas, fL4);
fMain->AddFrame(fF3, fL3);
fMain->AddFrame(fF5, fL4);
initHist();
drawGridlines();
drawSource();
std::string strtmp = "Total counts: " + std::to_string(counts);
countsText = new TText(0.7, 0.92, strtmp.c_str());
countsText->SetTextSizePixels(2);
countsText->SetNDC(kTRUE);
countsText->Draw("SAME");
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
// end embeded canvas
return 0;
}
void MyMainFrame::initHist() {
histo = new TH2D("ROI", " ; Azimuth; Elevation",
config->phiBins, -180, 180,
config->thetaBins, -90, 90);
// canvas = new TCanvas("THCanvas","Sinocanvas", 1000, 500);
// canvas->Connect("Closed()", "Image", this, "closed()");
// init image
for (int i = 0; i < config->phiBins; i++)
{
for (int j = 0; j < config->thetaBins; j++)
{
histo->SetBinContent(i+1, j+1, 0);
}
}
gStyle->SetOptStat(0);
histo->GetZaxis()->SetLabelSize(0.02);
histo->Draw("z aitoff");
}
void MyMainFrame::drawGridlines() {
float conv=M_PI/180;
float la, lo, x, y, z;
for (int j=0;j<Nl;++j) {
// latitudes[j]=new TGraph();
latitudes.push_back(new TGraph());
la = -90+180/(Nl-1)*j;
for (int i=0;i<M+1;++i) {
lo = -180+360/M*i;
z = sqrt(1+cos(la*conv)*cos(lo*conv/2));
x = 180*cos(la*conv)*sin(lo*conv/2)/z;
y = 90*sin(la*conv)/z;
latitudes[j]->SetPoint(i,x,y);
}
}
for (int j=0;j<NL;++j) {
// longitudes[j]=new TGraph();
longitudes.push_back(new TGraph());
lo = -180+360/(NL-1)*j;
for (int i=0;i<M+1;++i) {
la = -90+180/M*i;
z = sqrt(1+cos(la*conv)*cos(lo*conv/2));
x = 180*cos(la*conv)*sin(lo*conv/2)/z;
y = 90*sin(la*conv)/z;
longitudes[j]->SetPoint(i,x,y);
}
}
// Draw the grid
TPad *pad2 = new TPad("pad2","",0,0,1,1);
pad2->SetFillStyle(4000);
pad2->SetFillColor(0);
pad2->SetBorderSize(0);
pad2->SetFrameBorderMode(0);
pad2->SetFrameLineColor(0);
pad2->SetFrameBorderMode(0);
pad2->Draw();
pad2->cd();
Double_t ymin = -89.5;
Double_t ymax = 89.5;
Double_t dy = (ymax-ymin)/0.8; //10 per cent margins top and bottom
Double_t xmin = -180;
Double_t xmax = 180;
Double_t dx = (xmax-xmin)/0.8; //10 per cent margins left and right
pad2->Range(xmin-0.1*dx,ymin-0.1*dy,xmax+0.1*dx,ymax+0.1*dy);
for (int j=0;j<Nl;++j) latitudes[j]->Draw("l");
for (int j=0;j<NL;++j) longitudes[j]->Draw("l");
}
void MyMainFrame::drawSource() {
// draw source marker
float conv=M_PI/180;
double sigma = 1;
double Ax, Ay;
aitoff2xy(config->truePhi / conv, config->trueTheta / conv, Ax, Ay);
sourceMarker = new TEllipse(Ax, Ay, sigma, sigma);
sourceMarker->SetFillColor(0);
sourceMarker->SetFillStyle(0);
sourceMarker->SetLineWidth(2);
sourceMarker->SetLineColor(kRed);
sourceMarker->Draw("SAME");
}
void MyMainFrame::Start()
{
// printf("\"Start updating\" Button Pressed!\n");
stop = false;
}
void MyMainFrame::Stop()
{
// printf("\"Stop updating\" Button Pressed!\n");
stop = true;
}
void MyMainFrame::Clear()
{
// printf("\"Clear image\" Button Pressed!\n");
counts = 0;
histo->Reset();
redraw();
}
void MyMainFrame::addCone(std::vector<Cone>::const_iterator first,
std::vector<Cone>::const_iterator last)
{
if (first == last)
{
return;
}
// project cones onto the spherical surface
// std::cout << "Projecting cones onto the designated spherical surface..."<<std::endl;
double alpha(0);
double beta(0);
double sgma2(0);
double sgmb2(0);
const double sgmE(0.06); // 6% energy resolution
const Vector3D sgmpos(0.8 / 2, 0.43 / 2, 1.0 / 2); // mm
// const Vector3D sgmpos(0, 0, 0); // mm
const int order(3);
Vector3D ray;
sgma2 = std::pow(0.511*k->E0/std::pow((k->E0-k->Edpst),2)*sgmE, 2);
#pragma omp parallel for private(ray, beta, sgmb2)
for (int i = 0; i < config->thetaBins; i++)
{
for (int j = 0; j < config->phiBins; j++)
{
ray = k->apex - config->xbs[i][j];
beta = getCosAngle(ray, k->axis);
sgmb2 = std::pow((ray/(ray*ray) + k->axis/(k->axis*k->axis)-
(ray+k->axis)/(ray*k->axis))*sgmpos * beta, 2);
sgmb2+= std::pow((ray/(ray*k->axis)-k->axis/(k->axis*k->axis))*sgmpos * beta, 2);
histo->SetBinContent(j+1, i+1,
histo->GetBinContent(j+1, i+1) +
std::exp(-std::pow((std::pow(beta, order)-std::pow(alpha, order)), 2)/
(2*order*order*(std::pow(alpha,2*order-2)*sgma2 +std::pow(beta, 2*order-2)*sgmb2))
)
);
std::string strtmp = "Total counts: " + std::to_string(counts);
countsText->SetText(0.7, 0.92, strtmp.c_str());
canvas->Modified();
canvas->Update();
}
int aitoff2xy(const double& l, const double& b, double &Al, double &Ab)
{
Double_t x, y;
Double_t alpha2 = (l/2)*TMath::DegToRad();
Double_t delta = b*TMath::DegToRad();
Double_t r2 = TMath::Sqrt(2.);
Double_t f = 2*r2/TMath::Pi();
Double_t cdec = TMath::Cos(delta);
Double_t denom = TMath::Sqrt(1. + cdec*TMath::Cos(alpha2));
x = cdec*TMath::Sin(alpha2)*2.*r2/denom;
y = TMath::Sin(delta)*r2/denom;
x *= TMath::RadToDeg()/f;
y *= TMath::RadToDeg()/f;
// x *= -1.; // for a skymap swap left<->right
Al = x;
Ab = y;
return 0;
}