From ae69801404023b6305de124b7e6ef956c62fdbea Mon Sep 17 00:00:00 2001 From: mingf2 <fm140905@gmail.com> Date: Wed, 19 May 2021 01:42:54 -0500 Subject: [PATCH] normalize image --- .vscode/c_cpp_properties.json | 4 +-- Headers/gui.h | 3 ++ Headers/setup.h | 6 ++-- Sources/gui.cpp | 58 +++++++++++++++++++++++++++++++++++ Sources/main.cpp | 3 +- Sources/setup.cpp | 7 +++++ 6 files changed, 76 insertions(+), 5 deletions(-) diff --git a/.vscode/c_cpp_properties.json b/.vscode/c_cpp_properties.json index fa67db6..4a3ed5d 100644 --- a/.vscode/c_cpp_properties.json +++ b/.vscode/c_cpp_properties.json @@ -4,8 +4,8 @@ "name": "Linux", "includePath": [ "${workspaceFolder}/**", - "/home/ming/root/include", - "/home/ming/root/lib" + "/home/mingf2/root/include", + "/home/mingf2/root/lib" ], "defines": [], "compilerPath": "/usr/bin/gcc", diff --git a/Headers/gui.h b/Headers/gui.h index 4293573..976052f 100644 --- a/Headers/gui.h +++ b/Headers/gui.h @@ -19,6 +19,7 @@ #include <TGComboBox.h> #include <TGFileDialog.h> #include <TText.h> +#include <TPaletteAxis.h> #include "setup.h" @@ -63,6 +64,8 @@ public: MyMainFrame(const Setup* config_, const TGWindow *p,UInt_t w,UInt_t h); virtual ~MyMainFrame(); void addCone(std::vector<Cone>::const_iterator first, std::vector<Cone>::const_iterator last); + void addConeNormalized(std::vector<Cone>::const_iterator first, + std::vector<Cone>::const_iterator last); // slots void Start(); void Stop(); diff --git a/Headers/setup.h b/Headers/setup.h index ee3310b..6cb83eb 100644 --- a/Headers/setup.h +++ b/Headers/setup.h @@ -63,7 +63,7 @@ public: } Setup(const Vector3D v): R(std::sqrt(v.X*v.X + v.Y*v.Y + v.Z * v.Z)), - trueSource(v) + trueSource(v) { getPixelPos(); pixelateSphere(); @@ -78,10 +78,12 @@ public: const double dtheta = 180 / double(thetaBins) * M_PI / 180; const double thetaMin = -90 * M_PI / 180 + dtheta / 2; const double thetaMax = 90 * M_PI / 180 - dtheta / 2; + std::vector<double> thetaBinCenters; // Azimuthal angle, rad const double dphi = 2 * M_PI /double(phiBins); const double phiMin = -M_PI + dphi / 2; const double phiMax = M_PI - dphi / 2; + std::vector<double> phiBinCenters; // pixelate the projection sphere std::vector<std::vector<Vector3D>> xbs; int pixelateSphere(); @@ -98,7 +100,7 @@ public: const int order=3; // refresh every 10 events - const int chuckSize = 100; + const int chuckSize = 10; // queue capacity ulong capacity=1024; diff --git a/Sources/gui.cpp b/Sources/gui.cpp index ada7dde..138b5be 100644 --- a/Sources/gui.cpp +++ b/Sources/gui.cpp @@ -233,6 +233,64 @@ void MyMainFrame::addCone(std::vector<Cone>::const_iterator first, } +void MyMainFrame::addConeNormalized(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); + double summ(0); + Vector3D ray; + std::vector<std::vector<double>> probDist(config->thetaBins, std::vector<double>(config->phiBins, 0)); + for (auto k = first; k < last; k++) + { + counts++; + // update the image + alpha = k->cosHalfAngle; + sgma2 = std::pow(0.511*k->E0/std::pow((k->E0-k->Edpst),2)*config->sgmE, 2); + summ = 0; + #pragma omp parallel for private(ray, beta, sgmb2) shared(probDist) reduction(+:summ) + 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))*config->sgmpos * beta, 2); + sgmb2+= std::pow((ray/(ray*k->axis)-k->axis/(k->axis*k->axis))*config->sgmpos * beta, 2); + probDist[i][j] = std::exp(-std::pow((std::pow(beta, config->order)-std::pow(alpha, config->order)), 2)/ + (2*std::pow(config->order, 2)* + (std::pow(alpha,2*config->order-2)*sgma2 +std::pow(beta, 2*config->order-2)*sgmb2))); + // use proability density + // uncomment to use proability counts (integral) + // probDist[i][j] = probDist[i][j] * config->dtheta * config->dphi * std::cos(config->thetaBinCenters[i]); + // summ += probDist[i][j]; + summ += probDist[i][j] * config->dtheta * config->dphi * std::cos(config->thetaBinCenters[i]); + } + } + // #pragma omp parallel for + for (int i = 0; i < config->thetaBins; i++) + { + for (int j = 0; j < config->phiBins; j++) + { + histo->SetBinContent(j+1, i+1, (histo->GetBinContent(j+1, i+1)*(counts - 1) + probDist[i][j] / summ) / counts); + } + } + } + // counts += (last - first); + // redraw + redraw(); +} + void MyMainFrame::redraw() { std::string strtmp = "Total counts: " + std::to_string(counts); countsText->SetText(0.7, 0.92, strtmp.c_str()); diff --git a/Sources/main.cpp b/Sources/main.cpp index 93731ed..43a228f 100644 --- a/Sources/main.cpp +++ b/Sources/main.cpp @@ -38,7 +38,8 @@ int main(int argc, char** argv) { i++; } // update image - myMainFrame->addCone(cones.cbegin(), cones.cbegin() + i); + // myMainFrame->addCone(cones.cbegin(), cones.cbegin() + i); + myMainFrame->addConeNormalized(cones.cbegin(), cones.cbegin() + i); if (!unfinished) { // stop updating after this final run diff --git a/Sources/setup.cpp b/Sources/setup.cpp index ea628cf..18f6645 100644 --- a/Sources/setup.cpp +++ b/Sources/setup.cpp @@ -21,12 +21,19 @@ int Setup::getPixelPos() int Setup::pixelateSphere() { xbs.resize(thetaBins, std::vector<Vector3D>(phiBins, Vector3D(0, 0, 0))); + thetaBinCenters.resize(thetaBins, 0); + phiBinCenters.resize(phiBins, 0); double phi = phiMin; double theta = thetaMin; for (int i = 0; i < thetaBins; i++) { + thetaBinCenters[i] = theta; + phi = phiMin; for (int j = 0; j < phiBins; j++) { + if (i==0) + phiBinCenters[j] = phi; + xbs[i][j].X = R * std::cos(theta) * std::cos(phi); xbs[i][j].Y = R * std::cos(theta) * std::sin(phi); xbs[i][j].Z = R * std::sin(theta); -- GitLab