diff --git a/Postprocessing/.gitignore b/Postprocessing/.gitignore new file mode 100644 index 0000000000000000000000000000000000000000..48979898029c094e25cdd6bcd77ec4fd3ae6a9b4 --- /dev/null +++ b/Postprocessing/.gitignore @@ -0,0 +1 @@ +polimi \ No newline at end of file diff --git a/Postprocessing/CMakeLists.txt b/Postprocessing/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..bccf85a05f60236874d581f8a15e4e1c666daf4b --- /dev/null +++ b/Postprocessing/CMakeLists.txt @@ -0,0 +1,13 @@ +cmake_minimum_required(VERSION 3.16.3) + +project(Postprocessing CXX) + +set(CMAKE_CXX_STANDARD 11) +set(CMAKE_POSITION_INDEPENDENT_CODE ON) + +include_directories(${CMAKE_CURRENT_SOURCE_DIR}/Header) + +enable_testing() +add_subdirectory(googletest) +add_subdirectory(Sources) +add_subdirectory(Test) \ No newline at end of file diff --git a/Postprocessing/Header/reader.h b/Postprocessing/Header/reader.h new file mode 100644 index 0000000000000000000000000000000000000000..545f87b2e43e3e2c2b8bf488b0f46b3109884746 --- /dev/null +++ b/Postprocessing/Header/reader.h @@ -0,0 +1,115 @@ +#pragma once + +#include <sstream> +#include <string> +#include <vector> +#include <cmath> + +class Vector3D +{ +public: + double X=0; + double Y=0; + double Z=0; + Vector3D():Vector3D(0, 0, 0) {} + Vector3D(const double x, const double y, const double z): X(x), Y(y), Z(z) {} + + Vector3D operator+(const Vector3D& r) const { + return Vector3D(this->X + r.X, this->Y + r.Y, this->Z + r.Z); + } + Vector3D operator-(const Vector3D& r) const { + return Vector3D(this->X - r.X, this->Y - r.Y, this->Z - r.Z); + } + + double operator*(const Vector3D& r) const { + return this->X * r.X + this->Y * r.Y + this->Z * r.Z; + } + Vector3D operator*(const double& d) const { + return Vector3D(this->X * d, this->Y * d, this->Z * d); + } + Vector3D operator/(const double& d) const { + return Vector3D(this->X / d, this->Y / d, this->Z / d); + } +}; + +inline Vector3D operator*(const double& d, const Vector3D& r){ + return r * d; +} + +inline double getCosAngle(const Vector3D& l, const Vector3D& r) { + return l*r / std::sqrt((l*l) * (r*r)); +} + +class Cone +{ +public: + // interaction site 1 is the apex + Vector3D apex; + // site 2 - site 1 is the axis + Vector3D axis; + // half angle = abs (scattering angle) + // for neutron, 0 < half angle < pi /2 + // for gamma, 0 < half angle < pi + double cosHalfAngle; + // initial energy + double E0; + // energy deposited in compoton scattering + double Edpst; + Cone(): Cone(Vector3D(0,0,0), Vector3D(0,0,0), 0, 0, 0) {} + Cone(const Vector3D point, const Vector3D line, const double cosTheta, const double E0_, const double Edpst_): + apex(point), axis(line), cosHalfAngle(cosTheta), E0(E0_), Edpst(Edpst_) {} + Cone(const std::string& record) { + apex = Vector3D(std::stod(record.substr(0, 8)), + std::stod(record.substr(8, 8)), + std::stod(record.substr(16, 8))); + axis = Vector3D(std::stod(record.substr(24, 8)), + std::stod(record.substr(32, 8)), + std::stod(record.substr(40, 8))); + cosHalfAngle = std::stod(record.substr(48, 13)); + E0 = std::stod(record.substr(61, 13)); + Edpst = std::stod(record.substr(74, 13)); + } +}; + + +class Pulse +{ +private: + /* data */ +public: + ulong histNo=0; // history number + double height=0; // pulse amplitude, MeV + ushort reaction=0; // reaction type, 1 or 3 + int cellNo=0; // cell number + double ergIn=0; // energy of incoming photon, MeV + Pulse() {} + Pulse(std::string record) { + histNo = std::stoi(record.substr(56, 8)); + height = std::stod(record.substr(77, 9)); + reaction = std::stoi(record.substr(88, 4)); + cellNo = std::stoi(record.substr(137, 9)); + ergIn = std::stod(record.substr(161, 10)); + } +}; + +class Event +{ +private: + /* data */ +public: + ulong histNo=0; // history number + int cellNo=0; // cell number + double ergDpst=0; // energy of incoming photon, MeV + Vector3D pos=Vector3D(0,0,0); + Event() {} + Event(std::string record) { + histNo = std::stoi(record.substr(0, 10)); + cellNo = std::stoi(record.substr(35, 6)); + ergDpst = std::stod(record.substr(44, 8)); + pos = Vector3D(stod(record.substr(72,6)), + stod(record.substr(80,6)), + stod(record.substr(88,6))); + } +}; + +Vector3D getPos(const std::vector<Event>& events, const int& cellNo); \ No newline at end of file diff --git a/Postprocessing/Sources/CMakeLists.txt b/Postprocessing/Sources/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..2e130dacbfb14508c9eaa5a40beae2616bc4025c --- /dev/null +++ b/Postprocessing/Sources/CMakeLists.txt @@ -0,0 +1,4 @@ +add_library(reader SHARED reader.cpp) + +add_executable(main main.cpp) +target_link_libraries(main PUBLIC reader) \ No newline at end of file diff --git a/Postprocessing/Sources/main.cpp b/Postprocessing/Sources/main.cpp new file mode 100644 index 0000000000000000000000000000000000000000..904c3af0969267c216d0c75765124b916760b333 --- /dev/null +++ b/Postprocessing/Sources/main.cpp @@ -0,0 +1,175 @@ +#include <iostream> +#include <fstream> +#include <vector> +#include <string> +#include <algorithm> +#include <iomanip> +#include <chrono> +#include <set> + +#include "reader.h" + +int main(int argc, char** argv) +{ + // std::string dumnFile("/media/ming/DATA/projects/Imager/polimi/dumn1"); + // std::string pulseFile("/media/ming/DATA/projects/Imager/polimi/imager_All_pulses.o"); + // std::string outpath("/media/ming/DATA/projects/Imager/polimi/cones.txt"); + // test run + std::string dumnFile("/media/ming/DATA/projects/Imager/Test/collision.txt"); + std::string pulseFile("/media/ming/DATA/projects/Imager/Test/pulses.txt"); + std::string outpath("/media/ming/DATA/projects/Imager/Test/cones.txt"); + std::set<int> channels = {13, 15, 17, 19, 21, + 53, 55, 57, 59, 61, + 73, 75, 77, 79, 81, + 93, 95, 97, 99, 101, + 303, 305, 307, 309, 311, + 503, 505, 507, 509, 511}; + + std::ifstream fdumn; + fdumn.open(dumnFile, std::ios::in); + if (!fdumn.good()) + { + throw std::invalid_argument("Cannot open collision file"); + } + + std::ifstream fpulse; + fpulse.open(pulseFile, std::ios::in); + if (!fpulse.good()) + { + throw std::invalid_argument("Cannot open pulse file"); + } + + std::ofstream outfile; + outfile.open(outpath, std::ios::out); + if (!outfile.good()) + { + throw std::invalid_argument("Cannot create file: " + outpath); + } + + auto startTime = std::chrono::high_resolution_clock::now(); + // reconstruct events + std::vector<Cone> cones; + std::vector<ulong> histryNumbers; + double cosHalfangle; + Vector3D apex; + Vector3D axis; + std::string pulseLine; + std::string dumnLine; + std::string prevDumnLine(""); + std::vector<Pulse> pulses; // pulses with same history number + u_long histNo(0); + while (std::getline(fpulse, pulseLine)) + { + // parse line + Pulse newPulse(pulseLine); + // if line.histNo == histNo, add this pulse, go to next line + if (newPulse.histNo == histNo || histNo == 0 || pulses.size() == 0) + { + // // pick up channels of interest + // if (channels.find(newPulse.cellNo) != channels.end()) + { + // update histNo; + histNo = newPulse.histNo; + pulses.push_back(std::move(newPulse)); + continue; + } + + } + // else, process vector of pulses; + if (pulses.size() > 1) + { + // sort the pulses based on incoming photon energy + std::sort(pulses.begin(), pulses.end(), + [](const Pulse & a, const Pulse & b) -> bool + { + return a.ergIn > b.ergIn; + }); + // remove anomalies + if (std::abs(pulses[0].ergIn - 0.661) < 0.001 && pulses[0].height < pulses[0].ergIn) + { + // compute scattering angle + cosHalfangle = 1 - pulses[0].height * 0.511 / (pulses[0].ergIn * (pulses[0].ergIn - pulses[0].height)); + if (std::abs(cosHalfangle) < 1) + { + // discard cost >= 1, i.e., multiple compton scattering + // get position from collision file + std::vector<Event> interactions; + if (prevDumnLine.length() > 0) + { + u_long histNotmp = stoi(prevDumnLine.substr(0, 10)); + if (histNotmp == histNo) + { + interactions.push_back(Event(prevDumnLine)); + } + } + + while (std::getline(fdumn, dumnLine)) + { + u_long histNotmp = stoi(dumnLine.substr(0, 10)); + if (histNotmp < histNo) + continue; + else if (histNotmp > histNo) + { + // overshoot + prevDumnLine = dumnLine; + // process vector interactions + // find the interaction that deposits most energy in given cell + // apex = site 1 + // axis = site 2 -site 1 + apex = getPos(interactions, pulses[0].cellNo); + axis = getPos(interactions, pulses[1].cellNo) - apex; + break; + } + else + { + interactions.push_back(Event(dumnLine)); + } + } + if (fdumn.eof()) + { + break; + } + if (axis * axis > 0) + { + // axis could be a zero vector if interactions occur at the edge + cones.push_back(Cone(apex, axis, cosHalfangle, pulses[0].ergIn, pulses[0].height)); + histryNumbers.push_back(histNo); + } + } + } + } + // empty vector pulses; + pulses.clear(); + // update histNo; + histNo = newPulse.histNo; + // add current line to vector; + pulses.push_back(Pulse(std::move(newPulse))); + //go to next line + } + + fdumn.close(); + fpulse.close(); + + // save cones to text file + outfile << " apex.x apex.y apex.z axis.x axis.y axis.z cosHalfAngle E0 Edpst history_number\n"; + for (int i = 0; i < cones.size(); i++) + { + outfile << std::fixed << std::setprecision(2) + << std::setw(8) << cones[i].apex.X // << " " + << std::setw(8) << cones[i].apex.Y // << " " + << std::setw(8) << cones[i].apex.Z // << " " + << std::setw(8) << cones[i].axis.X // << " " + << std::setw(8) << cones[i].axis.Y // << " " + << std::setw(8) << cones[i].axis.Z // << " " + << std::fixed << std::setprecision(8) + << std::setw(13) << cones[i].cosHalfAngle + << std::setw(13) << cones[i].E0 + << std::setw(13) << cones[i].Edpst + << std::setw(13) << histryNumbers[i] << '\n'; + } + outfile.close(); + + auto endTime = std::chrono::high_resolution_clock::now(); + std::cout << std::chrono::duration_cast<std::chrono::milliseconds>(endTime - startTime).count() << "ms" << std::endl; + return 0; +} \ No newline at end of file diff --git a/Postprocessing/Sources/reader.cpp b/Postprocessing/Sources/reader.cpp new file mode 100644 index 0000000000000000000000000000000000000000..c50e41d4afb6f24e77a7adfa1fb2dc7e443bb644 --- /dev/null +++ b/Postprocessing/Sources/reader.cpp @@ -0,0 +1,18 @@ +#include "reader.h" + +Vector3D getPos(const std::vector<Event>& events, const int& cellNo) +{ + // the interaction that deposits most energy in the target cell + Vector3D pos; + double maxDpstErg(0); + for (int i = 0; i < events.size(); i++) + { + if (events[i].cellNo == cellNo && events[i].ergDpst > maxDpstErg) + { + maxDpstErg = events[i].ergDpst; + pos = events[i].pos; + } + } + return pos; + +} \ No newline at end of file diff --git a/Postprocessing/Test/CMakeLists.txt b/Postprocessing/Test/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..2c5e7d9949dc94e4385a629edeb1d4e9b1d96aa3 --- /dev/null +++ b/Postprocessing/Test/CMakeLists.txt @@ -0,0 +1,8 @@ +add_executable(readerTest test.cpp) +target_link_libraries(readerTest PUBLIC gtest_main reader) + +add_test( + NAME readerTest + COMMAND readerTest +) + \ No newline at end of file diff --git a/Postprocessing/Test/collision.txt b/Postprocessing/Test/collision.txt new file mode 100644 index 0000000000000000000000000000000000000000..8f8cd49d4c07b72c47efec704c890b33e1bbc6e0 --- /dev/null +++ b/Postprocessing/Test/collision.txt @@ -0,0 +1,200 @@ + 21 1 2 1 55 35 0.228927 0.04 0.84 2.26 -0.16 1.000E+00 0 0 0 6.610E-01 + 21 1 2 1 55 35 0.199572 0.04 1.11 2.57 -0.04 1.000E+00 0 1 0 4.220E-01 + 21 1 2 3 55 35 0.186618 0.04 1.42 2.27 -0.11 1.000E+00 0 2 0 2.224E-01 + 21 1 2 3 53 35 0.035813 0.04 1.42 2.27 -0.11 1.000E+00 0 2 0 3.581E-02 + 21 2 2 3 53 35 0.010071 0.04 0.84 2.26 -0.16 1.000E+00 0 1 0 1.007E-02 + 48 1 2 1 55 56 0.429851 0.04 1.76 2.62 -0.05 1.000E+00 0 0 0 6.610E-01 + 48 1 2 1 55 55 0.073002 0.04 1.95 2.45 0.11 1.000E+00 0 1 0 2.096E-01 + 48 1 2 3 55 55 0.105975 0.04 1.96 2.46 0.07 1.000E+00 0 2 0 1.366E-01 + 48 1 2 3 55 55 0.030622 0.04 1.97 2.46 0.08 1.000E+00 0 2 0 3.062E-02 + 48 3 2 3 55 56 0.013385 0.04 1.76 2.62 -0.05 1.000E+00 0 1 0 1.338E-02 + 48 2 2 3 55 56 0.008165 0.04 1.76 2.62 -0.05 1.000E+00 0 1 0 8.165E-03 + 51 1 2 1 53 74 0.050400 0.04 2.86 1.86 0.28 1.000E+00 0 1 0 6.281E-01 + 54 1 2 1 55 307 0.358192 0.04 5.11 3.38 0.52 1.000E+00 0 0 0 6.610E-01 + 54 1 2 1 53 98 0.094290 0.05 3.31 3.82 -0.28 1.000E+00 0 1 0 3.028E-01 + 54 1 2 3 53 77 0.179909 0.05 2.91 3.41 0.01 1.000E+00 0 2 0 2.085E-01 + 54 1 2 3 55 77 0.028609 0.05 2.90 3.40 0.00 1.000E+00 0 2 0 2.861E-02 + 59 1 2 1 53 604 0.108711 0.04 7.30 2.09 1.47 1.000E+00 0 0 0 6.610E-01 + 64 1 2 1 55 611 0.099624 0.05 7.21 5.02 -0.33 1.000E+00 0 1 0 6.588E-01 + 67 1 2 1 55 53 0.186788 0.04 2.16 1.40 -0.55 1.000E+00 0 0 0 6.610E-01 + 67 1 2 1 55 53 0.304438 0.04 2.18 1.67 -0.86 1.000E+00 0 1 0 4.742E-01 + 69 1 2 1 53 403 0.259676 0.04 6.30 1.52 -0.90 1.000E+00 0 0 0 6.610E-01 + 69 1 2 1 55 504 0.060774 0.04 7.05 1.88 -0.35 1.000E+00 0 1 0 3.939E-01 + 69 1 2 3 55 504 0.333151 0.05 7.18 2.10 0.37 1.000E+00 0 2 0 3.332E-01 + 69 2 2 3 55 403 0.007398 0.04 6.30 1.52 -0.90 1.000E+00 0 1 0 7.398E-03 + 86 1 2 1 53 404 0.054758 0.04 5.82 1.89 -0.02 1.000E+00 0 1 0 5.012E-01 + 86 1 2 1 53 404 0.110038 0.04 5.75 1.99 -0.07 1.000E+00 0 2 0 4.464E-01 + 86 1 2 3 55 304 0.336370 0.05 4.84 1.99 -0.82 1.000E+00 0 3 0 3.364E-01 + 95 1 2 1 53 209 0.023207 0.05 4.51 4.26 0.82 1.000E+00 0 0 0 6.610E-01 + 95 1 2 3 53 210 0.609487 0.05 4.53 4.45 0.78 1.000E+00 0 1 0 6.378E-01 + 95 1 2 3 53 210 0.024306 0.05 4.52 4.47 0.77 1.000E+00 0 1 0 2.831E-02 + 95 1 2 3 55 210 0.004000 0.05 4.52 4.47 0.78 1.000E+00 0 1 0 4.000E-03 + 97 1 2 1 55 39 0.009287 0.05 1.57 3.96 -0.64 1.000E+00 0 0 0 6.610E-01 + 97 1 2 1 53 41 0.171675 0.05 1.29 4.88 -0.69 1.000E+00 0 1 0 6.517E-01 + 97 1 2 1 55 42 0.304254 0.05 0.97 5.52 -1.45 1.000E+00 0 2 0 4.800E-01 + 97 1 2 3 55 41 0.144805 0.06 1.12 4.75 -1.22 1.000E+00 0 3 0 1.758E-01 + 97 1 2 3 53 41 0.030979 0.06 1.12 4.75 -1.22 1.000E+00 0 3 0 3.098E-02 + 100 1 2 1 53 80 0.457661 0.05 2.84 4.49 1.03 1.000E+00 0 0 0 6.610E-01 + 100 1 2 4 53 79 0.168368 0.05 2.66 4.21 0.88 1.000E+00 0 1 0 2.010E-01 + 100 1 2 3 55 79 0.028609 0.05 2.66 4.21 0.88 1.000E+00 0 1 0 2.861E-02 + 100 2 2 3 53 79 0.004000 0.05 2.66 4.21 0.88 1.000E+00 0 1 0 4.000E-03 + 100 2 2 3 55 80 0.002361 0.05 2.84 4.49 1.03 1.000E+00 0 1 0 2.361E-03 + 109 1 2 1 53 505 0.250427 0.04 6.49 2.16 -0.46 1.000E+00 0 0 0 6.610E-01 + 109 1 2 1 53 308 0.112826 0.05 5.12 3.65 0.13 1.000E+00 0 1 0 4.106E-01 + 109 1 2 1 55 58 0.058798 0.06 2.13 3.60 -0.78 1.000E+00 0 2 0 2.977E-01 + 109 1 2 3 55 59 0.207970 0.06 1.86 3.95 -0.80 1.000E+00 0 3 0 2.389E-01 + 109 1 2 3 53 59 0.030979 0.06 1.86 3.95 -0.80 1.000E+00 0 3 0 3.098E-02 + 122 1 2 1 55 98 0.241758 0.05 3.40 3.81 0.56 1.000E+00 0 0 0 6.610E-01 + 122 1 2 1 55 98 0.142120 0.05 3.33 3.86 0.52 1.000E+00 0 1 0 4.192E-01 + 122 1 2 1 53 77 0.140802 0.05 2.70 3.28 -0.18 1.000E+00 0 2 0 2.771E-01 + 122 1 2 3 55 77 0.105341 0.05 2.78 3.37 -0.15 1.000E+00 0 3 0 1.363E-01 + 122 1 2 3 53 77 0.030979 0.05 2.78 3.37 -0.14 1.000E+00 0 3 0 3.098E-02 + 135 1 2 1 55 42 0.256333 0.05 1.10 5.38 -0.11 1.000E+00 0 0 0 6.610E-01 + 146 1 2 1 55 37 0.418347 0.04 1.38 3.11 1.63 1.000E+00 0 0 0 6.610E-01 + 146 1 2 1 55 37 0.098123 0.05 1.53 3.09 1.51 1.000E+00 0 1 0 2.427E-01 + 146 1 2 3 55 36 0.113908 0.05 1.40 2.98 1.49 1.000E+00 0 2 0 1.445E-01 + 146 1 2 3 55 36 0.026250 0.05 1.39 2.97 1.48 1.000E+00 0 2 0 3.062E-02 + 146 1 2 3 55 36 0.004371 0.05 1.39 2.97 1.48 1.000E+00 0 2 0 4.371E-03 + 147 1 2 1 53 607 0.385736 0.05 7.59 3.07 0.68 1.000E+00 0 0 0 6.610E-01 + 147 2 2 3 55 607 0.032232 0.05 7.59 3.09 0.69 1.000E+00 0 1 0 6.285E-02 + 147 2 2 3 55 607 0.030622 0.05 7.54 3.11 0.66 1.000E+00 0 1 0 3.062E-02 + 163 1 2 1 55 203 0.270781 0.04 4.54 1.39 -0.90 1.000E+00 0 0 0 6.610E-01 + 163 1 2 1 53 203 0.157595 0.04 4.33 1.61 -1.33 1.000E+00 0 1 0 3.867E-01 + 163 2 2 3 53 203 0.003471 0.04 4.54 1.39 -0.90 1.000E+00 0 1 0 3.471E-03 + 164 1 2 1 53 304 0.056465 0.04 5.25 1.97 1.37 1.000E+00 0 0 0 6.610E-01 + 164 1 2 1 53 306 0.133486 0.04 5.39 2.66 1.75 1.000E+00 0 1 0 6.045E-01 + 164 1 2 1 55 509 0.201357 0.05 6.57 4.15 1.64 1.000E+00 0 2 0 4.710E-01 + 164 1 2 3 53 509 0.145724 0.05 6.52 4.13 1.63 1.000E+00 0 3 0 1.743E-01 + 164 1 2 3 55 509 0.028609 0.05 6.52 4.13 1.64 1.000E+00 0 3 0 2.861E-02 + 164 2 2 3 53 509 0.095358 0.05 6.58 4.16 1.64 1.000E+00 0 3 0 9.536E-02 + 176 1 2 1 53 37 0.461906 0.04 1.54 3.17 -0.58 1.000E+00 0 0 0 6.610E-01 + 176 1 2 3 55 36 0.168115 0.05 1.40 2.83 -0.62 1.000E+00 0 1 0 1.991E-01 + 176 1 2 3 55 36 0.030979 0.05 1.42 2.82 -0.65 1.000E+00 0 1 0 3.098E-02 + 178 1 2 1 53 204 0.235591 0.04 4.15 2.07 0.29 1.000E+00 0 0 0 6.610E-01 + 178 1 2 1 55 97 0.259918 0.05 3.27 3.09 -0.69 1.000E+00 0 1 0 4.254E-01 + 178 1 2 3 55 97 0.130518 0.05 3.45 3.03 -0.45 1.000E+00 0 2 0 1.655E-01 + 178 1 2 3 53 97 0.034973 0.05 3.45 3.03 -0.45 1.000E+00 0 2 0 3.497E-02 + 189 1 2 1 53 53 0.242956 0.04 2.03 1.43 0.25 1.000E+00 0 1 0 4.582E-01 + 189 2 2 3 55 53 0.002459 0.04 2.03 1.43 0.25 1.000E+00 0 2 0 2.459E-03 + 223 1 2 1 53 15 0.165790 0.04 0.38 2.45 -1.47 1.000E+00 0 0 0 6.610E-01 + 227 1 2 1 55 101 0.015550 0.05 3.71 5.07 1.35 1.000E+00 0 0 0 6.610E-01 + 230 1 2 1 53 95 0.445812 0.04 3.42 2.48 -1.00 1.000E+00 0 0 0 6.610E-01 + 230 1 2 3 55 74 0.180026 0.04 3.11 2.09 -1.27 1.000E+00 0 1 0 2.110E-01 + 230 1 2 3 53 74 0.030979 0.04 3.07 2.04 -1.26 1.000E+00 0 1 0 3.098E-02 + 230 2 2 3 53 95 0.004183 0.04 3.42 2.48 -1.00 1.000E+00 0 1 0 4.183E-03 + 231 1 2 1 55 305 0.280884 0.04 5.12 2.31 1.31 1.000E+00 0 0 0 6.610E-01 + 231 2 2 3 55 305 0.030801 0.04 5.13 2.33 1.30 1.000E+00 0 1 0 3.080E-02 + 233 1 2 1 53 73 0.094775 0.04 3.00 1.30 0.07 1.000E+00 0 1 0 3.057E-01 + 233 1 2 1 55 75 0.014118 0.04 2.46 2.32 1.32 1.000E+00 0 2 0 2.109E-01 + 235 1 2 1 53 93 0.103229 0.04 3.75 1.40 0.16 1.000E+00 0 0 0 6.610E-01 + 235 1 2 3 55 59 0.526792 0.05 2.35 4.06 -0.48 1.000E+00 0 1 0 5.578E-01 + 235 1 2 3 55 59 0.030979 0.05 2.35 4.07 -0.47 1.000E+00 0 1 0 3.098E-02 + 252 1 2 3 55 73 0.410970 0.04 2.57 1.33 -0.35 1.000E+00 0 2 0 4.153E-01 + 252 1 2 3 55 73 0.004371 0.04 2.57 1.33 -0.35 1.000E+00 0 2 0 4.371E-03 + 255 1 2 3 53 56 0.548844 0.04 2.00 2.82 -0.07 1.000E+00 0 0 0 6.610E-01 + 255 1 2 3 53 56 0.028609 0.04 2.01 2.82 -0.07 1.000E+00 0 0 0 2.861E-02 + 255 2 2 1 53 56 0.006062 0.04 1.99 2.80 -0.11 1.000E+00 0 0 0 8.355E-02 + 255 2 2 3 53 56 0.049179 0.04 1.98 2.81 -0.12 1.000E+00 0 1 0 7.749E-02 + 255 2 2 3 55 56 0.028306 0.04 1.96 2.81 -0.12 1.000E+00 0 1 0 2.831E-02 + 264 1 2 3 55 508 0.653891 0.05 6.56 3.58 -1.05 1.000E+00 0 0 0 6.610E-01 + 264 3 2 3 55 508 0.005380 0.05 6.56 3.58 -1.05 1.000E+00 0 0 0 5.380E-03 + 264 2 2 3 53 508 0.001728 0.05 6.56 3.58 -1.05 1.000E+00 0 0 0 1.728E-03 + 273 1 2 1 53 13 0.351409 0.04 0.55 1.44 -1.35 1.000E+00 0 1 0 6.297E-01 + 279 1 2 1 55 610 0.183755 0.05 7.38 4.66 1.31 1.000E+00 0 0 0 6.610E-01 + 293 1 2 1 53 307 0.301651 0.04 5.43 3.04 -1.60 1.000E+00 0 0 0 6.610E-01 + 295 1 2 3 55 308 0.630378 0.05 5.39 3.84 0.86 1.000E+00 0 0 0 6.610E-01 + 295 1 2 3 53 308 0.030622 0.05 5.41 3.81 0.83 1.000E+00 0 0 0 3.062E-02 + 301 1 2 1 55 507 0.153195 0.05 6.69 3.43 -0.32 1.000E+00 0 0 0 6.610E-01 + 301 1 2 1 53 611 0.076760 0.05 7.42 4.95 -1.71 1.000E+00 0 1 0 5.078E-01 + 318 1 2 1 55 33 0.021824 0.04 1.43 1.43 0.74 1.000E+00 0 0 0 6.610E-01 + 318 1 2 3 55 17 0.299911 0.05 0.66 3.33 1.12 1.000E+00 0 1 0 6.392E-01 + 318 1 2 3 53 17 0.030979 0.05 0.66 3.33 1.12 1.000E+00 0 1 0 3.098E-02 + 318 2 2 1 55 38 0.136722 0.05 0.84 3.62 1.22 1.000E+00 0 1 0 3.083E-01 + 318 2 2 1 53 37 0.021209 0.05 1.05 3.06 1.65 1.000E+00 0 2 0 1.716E-01 + 318 2 2 3 55 36 0.119376 0.05 0.93 2.91 1.72 1.000E+00 0 3 0 1.504E-01 + 320 1 2 1 55 53 0.106082 0.04 1.77 1.41 -0.38 1.000E+00 0 0 0 6.610E-01 + 320 1 2 1 55 55 0.371127 0.04 2.09 2.19 -0.37 1.000E+00 0 1 0 5.549E-01 + 320 1 2 3 53 33 0.183791 0.05 1.25 1.47 -0.68 1.000E+00 0 2 0 1.838E-01 + 329 1 2 1 53 403 0.440691 0.04 6.38 1.70 0.67 1.000E+00 0 0 0 6.610E-01 + 329 1 2 3 55 403 0.189330 0.04 6.19 1.59 0.82 1.000E+00 0 1 0 2.203E-01 + 329 1 2 3 55 403 0.030979 0.04 6.19 1.59 0.82 1.000E+00 0 1 0 3.098E-02 + 334 1 2 1 55 36 0.250791 0.04 1.55 2.84 1.63 1.000E+00 0 0 0 6.610E-01 + 334 2 2 3 55 36 0.019930 0.04 1.54 2.87 1.62 1.000E+00 0 1 0 1.993E-02 + 339 1 2 1 53 35 0.415561 0.04 1.05 2.42 1.33 1.000E+00 0 1 0 6.001E-01 + 339 1 2 3 53 34 0.156198 0.04 1.20 2.04 1.03 1.000E+00 0 2 0 1.845E-01 + 339 1 2 3 55 34 0.028306 0.05 1.27 2.01 1.03 1.000E+00 0 2 0 2.831E-02 + 342 1 2 1 53 74 0.181111 0.04 2.42 2.08 1.26 1.000E+00 0 0 0 6.610E-01 + 342 1 2 3 55 75 0.448910 0.04 2.46 2.45 1.69 1.000E+00 0 1 0 4.799E-01 + 342 1 2 3 55 75 0.030979 0.04 2.45 2.46 1.70 1.000E+00 0 1 0 3.098E-02 + 346 1 2 1 53 404 0.038808 0.04 5.95 1.91 -1.28 1.000E+00 0 0 0 6.610E-01 + 346 1 2 3 55 405 0.622192 0.04 5.93 2.19 -1.35 1.000E+00 0 1 0 6.222E-01 + 349 1 2 3 53 14 0.631001 0.04 0.08 1.81 1.52 1.000E+00 0 0 0 6.610E-01 + 349 1 2 3 55 14 0.028609 0.04 0.09 1.81 1.53 1.000E+00 0 0 0 2.861E-02 + 349 2 2 3 55 14 0.001389 0.04 0.08 1.81 1.52 1.000E+00 0 0 0 1.389E-03 + 351 1 2 1 55 97 0.421494 0.04 3.97 3.39 0.11 1.000E+00 0 0 0 6.610E-01 + 351 2 2 3 55 97 0.016296 0.04 3.97 3.39 0.11 1.000E+00 0 1 0 1.630E-02 + 360 1 2 1 53 403 0.204707 0.04 5.99 1.68 1.17 1.000E+00 0 1 0 5.914E-01 + 360 1 2 3 53 403 0.354371 0.04 5.92 1.71 1.35 1.000E+00 0 2 0 3.866E-01 + 360 1 2 3 55 404 0.032277 0.04 5.92 1.82 1.32 1.000E+00 0 2 0 3.228E-02 + 363 1 2 1 55 75 0.117440 0.04 2.82 2.35 -0.81 1.000E+00 0 0 0 6.610E-01 + 363 1 2 1 55 212 0.141444 0.05 4.52 5.31 -0.83 1.000E+00 0 1 0 5.436E-01 + 365 1 2 1 53 13 0.146701 0.04 0.45 1.67 0.05 1.000E+00 0 0 0 6.610E-01 + 365 2 2 3 55 13 0.031580 0.04 0.45 1.69 0.06 1.000E+00 0 1 0 6.256E-02 + 365 2 2 3 53 13 0.030979 0.04 0.45 1.69 0.07 1.000E+00 0 1 0 3.098E-02 + 371 1 2 1 53 307 0.382108 0.04 5.08 3.30 0.16 1.000E+00 0 0 0 6.610E-01 + 371 1 2 3 55 307 0.248099 0.05 5.31 3.20 0.19 1.000E+00 0 1 0 2.481E-01 + 371 2 2 3 53 307 0.026793 0.04 5.08 3.30 0.16 1.000E+00 0 1 0 3.079E-02 + 371 2 2 3 53 307 0.004000 0.04 5.08 3.30 0.16 1.000E+00 0 1 0 4.000E-03 + 388 1 2 1 55 13 0.116915 0.04 0.48 1.60 1.44 1.000E+00 0 0 0 6.610E-01 + 388 1 2 1 55 40 0.056205 0.05 1.14 4.44 0.62 1.000E+00 0 1 0 5.300E-01 + 388 2 2 3 53 13 0.010110 0.04 0.48 1.60 1.45 1.000E+00 0 1 0 1.411E-02 + 388 2 2 3 55 13 0.004000 0.04 0.48 1.60 1.45 1.000E+00 0 1 0 4.000E-03 + 403 1 2 1 55 310 0.401048 0.05 5.47 4.61 -0.03 1.000E+00 0 0 0 6.610E-01 + 403 1 2 1 53 310 0.059365 0.05 5.22 4.48 -0.79 1.000E+00 0 1 0 2.600E-01 + 403 1 2 1 55 410 0.003038 0.05 5.64 4.47 -1.05 1.000E+00 0 2 0 2.006E-01 + 403 1 2 1 53 410 0.030585 0.05 5.67 4.46 -1.06 1.000E+00 0 3 0 1.975E-01 + 403 1 2 3 53 410 0.134687 0.05 5.93 4.59 -0.91 1.000E+00 0 4 0 1.670E-01 + 403 1 2 3 55 410 0.032277 0.05 5.94 4.59 -0.90 1.000E+00 0 4 0 3.228E-02 + 411 1 2 1 55 80 0.173578 0.05 2.82 4.66 1.06 1.000E+00 0 0 0 6.610E-01 + 411 1 2 1 53 102 0.309355 0.05 3.35 5.50 1.70 1.000E+00 0 1 0 4.874E-01 + 411 1 2 3 55 101 0.143094 0.05 3.39 4.97 1.50 1.000E+00 0 2 0 1.781E-01 + 411 1 2 3 53 101 0.034973 0.05 3.38 4.97 1.51 1.000E+00 0 2 0 3.497E-02 + 425 1 2 1 55 53 0.407581 0.04 2.06 1.46 -1.00 1.000E+00 0 0 0 6.610E-01 + 426 1 2 1 53 603 0.448854 0.04 7.71 1.58 1.24 1.000E+00 0 0 0 6.610E-01 + 426 1 2 1 53 603 0.031061 0.04 7.72 1.50 1.30 1.000E+00 0 1 0 2.100E-01 + 426 2 2 3 53 603 0.002130 0.04 7.71 1.58 1.24 1.000E+00 0 1 0 2.130E-03 + 440 1 2 1 55 203 0.092391 0.04 4.09 1.57 0.12 1.000E+00 0 0 0 6.610E-01 + 440 1 2 1 55 77 0.028850 0.05 3.12 3.29 0.08 1.000E+00 0 1 0 5.686E-01 + 440 1 2 1 55 81 0.346691 0.05 2.57 4.90 -0.38 1.000E+00 0 2 0 5.398E-01 + 440 1 2 3 55 81 0.162090 0.05 2.82 4.79 -0.26 1.000E+00 0 3 0 1.931E-01 + 440 1 2 3 55 81 0.030979 0.05 2.83 4.79 -0.25 1.000E+00 0 3 0 3.098E-02 + 464 1 2 1 53 403 0.375619 0.04 5.67 1.50 -0.83 1.000E+00 0 0 0 6.610E-01 + 464 1 2 3 55 203 0.254402 0.04 4.77 1.66 -0.56 1.000E+00 0 1 0 2.854E-01 + 464 1 2 3 55 203 0.026607 0.04 4.73 1.70 -0.59 1.000E+00 0 1 0 3.098E-02 + 464 1 2 3 53 203 0.004371 0.04 4.73 1.70 -0.59 1.000E+00 0 1 0 4.371E-03 + 468 1 2 1 55 34 0.179628 0.04 1.21 1.85 -0.16 1.000E+00 0 1 0 5.636E-01 + 468 1 2 1 53 57 0.137205 0.05 1.90 3.06 -1.28 1.000E+00 0 2 0 3.840E-01 + 468 1 2 3 53 77 0.218467 0.05 2.53 3.27 -0.97 1.000E+00 0 3 0 2.468E-01 + 468 1 2 3 53 77 0.028306 0.05 2.52 3.26 -0.98 1.000E+00 0 3 0 2.831E-02 + 470 1 2 1 53 36 0.275545 0.04 0.98 2.82 0.90 1.000E+00 0 0 0 6.610E-01 + 472 1 2 1 53 38 0.159336 0.05 1.28 3.51 0.01 1.000E+00 0 1 0 5.193E-01 + 472 1 2 1 53 19 0.077187 0.05 0.01 4.15 0.30 1.000E+00 0 2 0 3.599E-01 + 474 1 2 1 55 505 0.230944 0.04 6.76 2.37 0.95 1.000E+00 0 0 0 6.610E-01 + 485 1 2 1 55 206 0.013824 0.04 4.44 2.82 0.00 1.000E+00 0 0 0 6.610E-01 + 492 1 2 3 55 96 0.593870 0.04 3.60 2.71 0.31 1.000E+00 0 0 0 6.610E-01 + 492 1 2 3 55 96 0.030622 0.04 3.60 2.71 0.31 1.000E+00 0 0 0 3.062E-02 + 492 2 2 3 53 96 0.007898 0.04 3.60 2.71 0.31 1.000E+00 0 0 0 3.651E-02 + 492 2 2 3 55 96 0.028609 0.04 3.60 2.70 0.31 1.000E+00 0 0 0 2.861E-02 + 500 1 2 1 53 17 0.195684 0.05 0.18 3.42 0.59 1.000E+00 0 0 0 6.610E-01 + 500 1 2 1 53 19 0.050466 0.05 0.16 4.25 -0.20 1.000E+00 0 1 0 4.653E-01 + 525 1 2 1 53 204 0.381636 0.04 4.35 1.76 1.57 1.000E+00 0 0 0 6.610E-01 + 551 1 2 1 53 405 0.288397 0.04 5.74 2.24 1.13 1.000E+00 0 0 0 6.610E-01 + 551 1 2 3 53 405 0.344296 0.04 5.84 2.27 1.08 1.000E+00 0 1 0 3.726E-01 + 551 1 2 3 55 405 0.028306 0.04 5.84 2.27 1.08 1.000E+00 0 1 0 2.831E-02 + 553 1 2 1 53 206 0.448776 0.04 4.48 2.75 -1.11 1.000E+00 0 0 0 6.610E-01 + 553 1 2 1 55 206 0.025213 0.04 4.48 2.60 -1.25 1.000E+00 0 1 0 2.080E-01 + 553 1 2 3 55 205 0.182825 0.04 4.47 2.54 -1.24 1.000E+00 0 2 0 1.828E-01 + 553 2 2 3 55 206 0.004186 0.04 4.48 2.75 -1.11 1.000E+00 0 1 0 4.186E-03 + 557 1 2 1 55 93 0.148389 0.04 3.48 1.59 -0.31 1.000E+00 0 0 0 6.610E-01 + 557 1 2 1 55 40 0.104767 0.05 0.91 4.46 -1.56 1.000E+00 0 1 0 4.868E-01 + 557 1 2 1 55 20 0.013422 0.05 0.67 4.50 -1.55 1.000E+00 0 2 0 3.821E-01 diff --git a/Postprocessing/Test/cones.txt b/Postprocessing/Test/cones.txt new file mode 100644 index 0000000000000000000000000000000000000000..507d76e471327471959116eabcf9783841d11875 --- /dev/null +++ b/Postprocessing/Test/cones.txt @@ -0,0 +1,25 @@ + apex.x apex.y apex.z axis.x axis.y axis.z cosHalfAngle E0 Edpst history_number + 1.76 2.62 -0.05 0.20 -0.16 0.12 -0.68048645 0.66100001 0.45273100 48 + 5.11 3.38 0.52 -1.80 0.44 -0.80 0.06025452 0.66100001 0.36266100 54 + 6.30 1.52 -0.90 0.88 0.58 1.27 0.51782044 0.66100001 0.25391000 69 + 4.51 4.26 0.82 0.02 0.19 -0.04 0.97239547 0.66100001 0.02278900 95 + 1.57 3.96 -0.64 -0.28 0.92 -0.05 0.99492261 0.66100001 0.00431300 97 + 2.84 4.49 1.03 -0.18 -0.28 -0.15 -0.61549360 0.66100001 0.44706400 100 + 6.49 2.16 -0.46 -1.37 1.49 0.59 0.51941873 0.66100001 0.25339100 109 + 3.40 3.81 0.56 -0.70 -0.53 -0.74 -0.00225151 0.66100001 0.37316500 122 + 5.25 1.97 1.37 0.14 0.69 0.38 0.92122555 0.66100001 0.06112600 164 + 1.54 3.17 -0.58 -0.14 -0.34 -0.04 -0.70780034 0.66100001 0.45502400 176 + 4.15 2.07 0.29 -0.88 1.02 -0.98 0.58314659 0.66100001 0.23156100 178 + 3.42 2.48 -1.00 -0.31 -0.39 -0.27 -0.52020828 0.66100001 0.43817500 230 + 3.75 1.40 0.16 -1.40 2.66 -0.64 0.84702101 0.66100001 0.10919400 235 + 6.69 3.43 -0.32 0.73 1.52 -1.39 0.78720127 0.66100001 0.14267600 301 + 1.43 1.43 0.74 -0.59 2.19 0.48 0.97727465 0.66100001 0.01887600 318 + 1.77 1.41 -0.38 0.32 0.78 0.01 0.85542739 0.66100001 0.10413900 320 + 2.42 2.08 1.26 0.04 0.37 0.43 0.69093463 0.66100001 0.18878600 342 + 5.95 1.91 -1.28 -0.02 0.28 -0.07 0.95282269 0.66100001 0.03801800 346 + 2.82 2.35 -0.81 1.70 2.96 -0.02 0.84884332 0.66100001 0.10810600 363 + 0.48 1.60 1.44 0.66 2.84 -0.82 0.79464811 0.66100001 0.13873100 388 + 5.47 4.61 -0.03 0.46 -0.02 -0.88 -0.78130296 0.66100001 0.46095100 403 + 2.82 4.66 1.06 0.53 0.84 0.64 0.72948559 0.66100001 0.17134200 411 + 4.09 1.57 0.12 -0.97 1.72 -0.04 0.85417659 0.66100001 0.10489700 440 + 5.67 1.50 -0.83 -0.90 0.16 0.27 -0.03122455 0.66100001 0.37778700 464 diff --git a/Postprocessing/Test/pulses.txt b/Postprocessing/Test/pulses.txt new file mode 100644 index 0000000000000000000000000000000000000000..c9b2d86a7e9cf4675bff3edb4d13cc47ab5f0da4 --- /dev/null +++ b/Postprocessing/Test/pulses.txt @@ -0,0 +1,100 @@ + 13. 0.40000000 2. 21. 0.737400 1. 1. -999. 0.10000000E+01 35 0.66100001 + 23. 0.40000000 2. 48. 0.215864 3. 1. -999. 0.10000000E+01 55 0.03062000 + 24. 0.40000000 2. 48. 0.452731 1. 1. -999. 0.10000000E+01 56 0.66100001 + 32. 0.40000000 2. 51. 0.046661 1. 1. -999. 0.10000000E+01 74 0.62809998 + 65. 0.40000000 2. 54. 0.362661 1. 1. -999. 0.10000000E+01 307 0.66100001 + 46. 0.50000000 2. 54. 0.083879 1. 1. -999. 0.10000000E+01 98 0.30280000 + 35. 0.50000000 2. 54. 0.207987 3. 1. -999. 0.10000000E+01 77 0.20850000 + 92. 0.40000000 2. 59. 0.103136 1. 1. -999. 0.10000000E+01 604 0.66100001 + 99. 0.50000000 2. 64. 0.102308 1. 1. -999. 0.10000000E+01 611 0.65880001 + 21. 0.40000000 2. 67. 0.503015 1. 1. -999. 0.10000000E+01 53 0.66100001 + 71. 0.40000000 2. 69. 0.253910 1. 1. -999. 0.10000000E+01 403 0.66100001 + 82. 0.40000000 2. 69. 0.390592 1. 1. -999. 0.10000000E+01 504 0.39390001 + 72. 0.40000000 2. 86. 0.172977 1. 1. -999. 0.10000000E+01 404 0.50120002 + 62. 0.50000000 2. 86. 0.341322 3. 1. -999. 0.10000000E+01 304 0.33640000 + 57. 0.50000000 2. 95. 0.022789 1. 1. -999. 0.10000000E+01 209 0.66100001 + 58. 0.50000000 2. 95. 0.631979 3. 1. -999. 0.10000000E+01 210 0.02831000 + 17. 0.50000000 2. 97. 0.004313 1. 1. -999. 0.10000000E+01 39 0.66100001 + 19. 0.50000000 2. 97. 0.351328 1. 1. -999. 0.10000000E+01 41 0.65170002 + 20. 0.50000000 2. 97. 0.295631 1. 1. -999. 0.10000000E+01 42 0.47999999 + 37. 0.50000000 2. 100. 0.207598 3. 2. -999. 0.10000000E+01 79 0.00400000 + 38. 0.50000000 2. 100. 0.447064 1. 1. -999. 0.10000000E+01 80 0.66100001 + 83. 0.40000000 2. 109. 0.253391 1. 1. -999. 0.10000000E+01 505 0.66100001 + 66. 0.50000000 2. 109. 0.117387 1. 1. -999. 0.10000000E+01 308 0.41060001 + 26. 0.60000000 2. 109. 0.061210 1. 1. -999. 0.10000000E+01 58 0.29769999 + 27. 0.60000000 2. 109. 0.240437 3. 1. -999. 0.10000000E+01 59 0.23890001 + 35. 0.50000000 2. 122. 0.286155 3. 1. -999. 0.10000000E+01 77 0.13630000 + 46. 0.50000000 2. 122. 0.373165 1. 1. -999. 0.10000000E+01 98 0.66100001 + 20. 0.50000000 2. 135. 0.253535 1. 1. -999. 0.10000000E+01 42 0.66100001 + 15. 0.40000000 2. 146. 0.523320 1. 1. -999. 0.10000000E+01 37 0.66100001 + 14. 0.50000000 2. 146. 0.137591 3. 1. -999. 0.10000000E+01 36 0.03062000 + 95. 0.50000000 2. 147. 0.456155 1. 1. -999. 0.10000000E+01 607 0.66100001 + 51. 0.40000000 2. 163. 0.425994 1. 1. -999. 0.10000000E+01 203 0.66100001 + 62. 0.40000000 2. 164. 0.061126 1. 1. -999. 0.10000000E+01 304 0.66100001 + 64. 0.40000000 2. 164. 0.139104 1. 1. -999. 0.10000000E+01 306 0.60450000 + 87. 0.50000000 2. 164. 0.489163 3. 1. -999. 0.10000000E+01 509 0.02861000 + 15. 0.40000000 2. 176. 0.455024 1. 1. -999. 0.10000000E+01 37 0.66100001 + 14. 0.50000000 2. 176. 0.190237 3. 1. -999. 0.10000000E+01 36 0.19910000 + 52. 0.40000000 2. 178. 0.231561 1. 1. -999. 0.10000000E+01 204 0.66100001 + 45. 0.50000000 2. 178. 0.452133 3. 1. -999. 0.10000000E+01 97 0.16550000 + 21. 0.40000000 2. 189. 0.253628 1. 1. -999. 0.10000000E+01 53 0.45820001 + 3. 0.40000000 2. 223. 0.167791 1. 1. -999. 0.10000000E+01 15 0.66100001 + 49. 0.50000000 2. 227. 0.020326 1. 1. -999. 0.10000000E+01 101 0.66100001 + 32. 0.40000000 2. 230. 0.207351 3. 1. -999. 0.10000000E+01 74 0.03098000 + 43. 0.40000000 2. 230. 0.438175 1. 1. -999. 0.10000000E+01 95 0.66100001 + 63. 0.40000000 2. 231. 0.311037 1. 1. -999. 0.10000000E+01 305 0.66100001 + 31. 0.40000000 2. 233. 0.091918 1. 1. -999. 0.10000000E+01 73 0.30570000 + 33. 0.40000000 2. 233. 0.013908 1. 1. -999. 0.10000000E+01 75 0.21089999 + 41. 0.40000000 2. 235. 0.109194 1. 1. -999. 0.10000000E+01 93 0.66100001 + 27. 0.50000000 2. 235. 0.553914 3. 1. -999. 0.10000000E+01 59 0.55779999 + 31. 0.40000000 2. 252. 0.407080 3. 1. -999. 0.10000000E+01 73 0.41530001 + 24. 0.40000000 2. 255. 0.683741 3. 1. -999. 0.10000000E+01 56 0.66100001 + 86. 0.50000000 2. 264. 0.701056 3. 1. -999. 0.10000000E+01 508 0.66100001 + 1. 0.40000000 2. 273. 0.363796 1. 1. -999. 0.10000000E+01 13 0.62970001 + 98. 0.50000000 2. 279. 0.183004 1. 1. -999. 0.10000000E+01 610 0.66100001 + 65. 0.40000000 2. 293. 0.289512 1. 1. -999. 0.10000000E+01 307 0.66100001 + 66. 0.50000000 2. 295. 0.671339 3. 1. -999. 0.10000000E+01 308 0.66100001 + 85. 0.50000000 2. 301. 0.142676 1. 1. -999. 0.10000000E+01 507 0.66100001 + 99. 0.50000000 2. 301. 0.069063 1. 1. -999. 0.10000000E+01 611 0.50779998 + 11. 0.40000000 2. 318. 0.018876 1. 1. -999. 0.10000000E+01 33 0.66100001 + 15. 0.50000000 2. 318. 0.021255 1. 2. -999. 0.10000000E+01 37 0.17160000 + 14. 0.50000000 2. 318. 0.126043 3. 2. -999. 0.10000000E+01 36 0.15040000 + 5. 0.50000000 2. 318. 0.324207 3. 1. -999. 0.10000000E+01 17 0.03098000 + 16. 0.50000000 2. 318. 0.140536 1. 2. -999. 0.10000000E+01 38 0.30829999 + 23. 0.40000000 2. 320. 0.363872 1. 1. -999. 0.10000000E+01 55 0.55489999 + 21. 0.40000000 2. 320. 0.104139 1. 1. -999. 0.10000000E+01 53 0.66100001 + 11. 0.50000000 2. 320. 0.183944 3. 1. -999. 0.10000000E+01 33 0.18380000 + 71. 0.40000000 2. 329. 0.669109 1. 1. -999. 0.10000000E+01 403 0.66100001 + 14. 0.40000000 2. 334. 0.273335 1. 1. -999. 0.10000000E+01 36 0.66100001 + 12. 0.40000000 2. 339. 0.186464 3. 1. -999. 0.10000000E+01 34 0.18449999 + 13. 0.40000000 2. 339. 0.418096 1. 1. -999. 0.10000000E+01 35 0.60009998 + 32. 0.40000000 2. 342. 0.188786 1. 1. -999. 0.10000000E+01 74 0.66100001 + 33. 0.40000000 2. 342. 0.451286 3. 1. -999. 0.10000000E+01 75 0.47990000 + 72. 0.40000000 2. 346. 0.038018 1. 1. -999. 0.10000000E+01 404 0.66100001 + 73. 0.40000000 2. 346. 0.605070 3. 1. -999. 0.10000000E+01 405 0.62220001 + 2. 0.40000000 2. 349. 0.649951 3. 1. -999. 0.10000000E+01 14 0.66100001 + 45. 0.40000000 2. 351. 0.435611 1. 1. -999. 0.10000000E+01 97 0.66100001 + 71. 0.40000000 2. 360. 0.570391 1. 1. -999. 0.10000000E+01 403 0.59140003 + 72. 0.40000000 2. 360. 0.040621 3. 1. -999. 0.10000000E+01 404 0.03228000 + 33. 0.40000000 2. 363. 0.108106 1. 1. -999. 0.10000000E+01 75 0.66100001 + 60. 0.50000000 2. 363. 0.139398 1. 1. -999. 0.10000000E+01 212 0.54360002 + 1. 0.40000000 2. 365. 0.212371 1. 1. -999. 0.10000000E+01 13 0.66100001 + 65. 0.40000000 2. 371. 0.671448 1. 1. -999. 0.10000000E+01 307 0.66100001 + 1. 0.40000000 2. 388. 0.138731 1. 1. -999. 0.10000000E+01 13 0.66100001 + 18. 0.50000000 2. 388. 0.060128 1. 1. -999. 0.10000000E+01 40 0.52999997 + 68. 0.50000000 2. 403. 0.460951 1. 1. -999. 0.10000000E+01 310 0.66100001 + 78. 0.50000000 2. 403. 0.207285 3. 1. -999. 0.10000000E+01 410 0.16700000 + 38. 0.50000000 2. 411. 0.171342 1. 1. -999. 0.10000000E+01 80 0.66100001 + 49. 0.50000000 2. 411. 0.164831 3. 1. -999. 0.10000000E+01 101 0.17810000 + 50. 0.50000000 2. 411. 0.306314 1. 1. -999. 0.10000000E+01 102 0.48740000 + 21. 0.40000000 2. 425. 0.400719 1. 1. -999. 0.10000000E+01 53 0.66100001 + 91. 0.40000000 2. 426. 0.494278 1. 1. -999. 0.10000000E+01 603 0.66100001 + 51. 0.40000000 2. 440. 0.104897 1. 1. -999. 0.10000000E+01 203 0.66100001 + 39. 0.50000000 2. 440. 0.561963 3. 1. -999. 0.10000000E+01 81 0.19310001 + 35. 0.50000000 2. 440. 0.028253 1. 1. -999. 0.10000000E+01 77 0.56860000 + 51. 0.40000000 2. 464. 0.280872 3. 1. -999. 0.10000000E+01 203 0.03098000 + 71. 0.40000000 2. 464. 0.377787 1. 1. -999. 0.10000000E+01 403 0.66100001 + 12. 0.40000000 2. 468. 0.185954 1. 1. -999. 0.10000000E+01 34 0.56360000 + 25. 0.50000000 2. 468. 0.149822 1. 1. -999. 0.10000000E+01 57 0.38400000 + 35. 0.50000000 2. 468. 0.248850 3. 1. -999. 0.10000000E+01 77 0.24680001 + 14. 0.40000000 2. 470. 0.270426 1. 1. -999. 0.10000000E+01 36 0.66100001 diff --git a/Postprocessing/Test/test.cpp b/Postprocessing/Test/test.cpp new file mode 100644 index 0000000000000000000000000000000000000000..f45b5f65988b20e2d399419892c135cea9e51b6c --- /dev/null +++ b/Postprocessing/Test/test.cpp @@ -0,0 +1,116 @@ +#include <gtest/gtest.h> +#include <string> + +#include "reader.h" + +TEST(PulseTest, constructor) +{ + std::string record(" 13. 0.40000000 2. 21. 0.737400 1. 1. -999. 0.10000000E+01 35 0.66100001"); + Pulse pulse(record); + EXPECT_EQ(pulse.histNo, 21); + EXPECT_DOUBLE_EQ(pulse.height, 0.737400); + EXPECT_EQ(pulse.reaction, 1); + EXPECT_EQ(pulse.cellNo, 35); + EXPECT_DOUBLE_EQ(pulse.ergIn, 0.66100001); + + record = " 22. 0.40000000 2. 89765585. 0.114611 3. 1. -999. 0.10000000E+01 54 0.11270000"; + pulse = Pulse(record); + EXPECT_EQ(pulse.histNo, 89765585); + EXPECT_DOUBLE_EQ(pulse.height, 0.114611); + EXPECT_EQ(pulse.reaction, 3); + EXPECT_EQ(pulse.cellNo, 54); + EXPECT_DOUBLE_EQ(pulse.ergIn, 0.11270000); +} + +TEST(EventTest, constructor) +{ + std::string record(" 21 1 2 1 55 35 0.228927 0.04 0.84 2.26 -0.16 1.000E+00 0 0 0 6.610E-01"); + Event event(record); + EXPECT_EQ(event.histNo, 21); + EXPECT_EQ(event.cellNo, 35); + EXPECT_DOUBLE_EQ(event.ergDpst, 0.228927); + Vector3D pos(0.84, 2.26, -0.16); + EXPECT_DOUBLE_EQ(event.pos.X, pos.X); + EXPECT_DOUBLE_EQ(event.pos.Y, pos.Y); + EXPECT_DOUBLE_EQ(event.pos.Z, pos.Z); + + record = " 99999990 1 2 1 55 604 0.202758 0.04 8.00 1.72 -0.13 1.000E+00 0 0 0 6.610E-01"; + event = Event(record); + EXPECT_EQ(event.histNo, 99999990); + EXPECT_EQ(event.cellNo, 604); + EXPECT_DOUBLE_EQ(event.ergDpst, 0.202758 ); + pos = Vector3D(8.00, 1.72, -0.13); + EXPECT_DOUBLE_EQ(event.pos.X, pos.X); + EXPECT_DOUBLE_EQ(event.pos.Y, pos.Y); + EXPECT_DOUBLE_EQ(event.pos.Z, pos.Z); +} + +TEST(ConeTest, constructor) +{ + std::string record(" 1.76 2.62 -0.05 0.20 -0.16 0.12 -0.68048645 0.66100001 0.45273100 48"); + Cone cone(record); + EXPECT_DOUBLE_EQ(cone.apex.X, 1.76); + EXPECT_DOUBLE_EQ(cone.apex.Y, 2.62); + EXPECT_DOUBLE_EQ(cone.apex.Z, -0.05); + EXPECT_DOUBLE_EQ(cone.axis.X, 0.20); + EXPECT_DOUBLE_EQ(cone.axis.Y, -0.16); + EXPECT_DOUBLE_EQ(cone.axis.Z, 0.12); + EXPECT_DOUBLE_EQ(cone.cosHalfAngle, -0.68048645); + EXPECT_DOUBLE_EQ(cone.E0, 0.66100001); + EXPECT_DOUBLE_EQ(cone.Edpst, 0.45273100); + + record = " 1.43 1.43 0.74 -0.59 2.19 0.48 0.97727465 0.66100001 0.01887600 318"; + cone = Cone(record); + EXPECT_DOUBLE_EQ(cone.apex.X, 1.43); + EXPECT_DOUBLE_EQ(cone.apex.Y, 1.43); + EXPECT_DOUBLE_EQ(cone.apex.Z, 0.74); + EXPECT_DOUBLE_EQ(cone.axis.X, -0.59); + EXPECT_DOUBLE_EQ(cone.axis.Y, 2.19); + EXPECT_DOUBLE_EQ(cone.axis.Z, 0.48); + EXPECT_DOUBLE_EQ(cone.cosHalfAngle, 0.97727465); + EXPECT_DOUBLE_EQ(cone.E0, 0.66100001); + EXPECT_DOUBLE_EQ(cone.Edpst, 0.01887600); +} + +TEST(getPosTest, getPosTest){ + std::vector<std::string> records = { + " 48 1 2 1 55 56 0.429851 0.04 1.76 2.62 -0.05 1.000E+00 0 0 0 6.610E-01", + " 48 1 2 1 55 55 0.073002 0.04 1.95 2.45 0.11 1.000E+00 0 1 0 2.096E-01", + " 48 1 2 3 55 55 0.105975 0.04 1.96 2.46 0.07 1.000E+00 0 2 0 1.366E-01", + " 48 1 2 3 55 55 0.030622 0.04 1.97 2.46 0.08 1.000E+00 0 2 0 3.062E-02", + " 48 3 2 3 55 56 0.013385 0.04 1.76 2.62 -0.05 1.000E+00 0 1 0 1.338E-02", + " 48 2 2 3 55 56 0.008165 0.04 1.76 2.62 -0.05 1.000E+00 0 1 0 8.165E-03"}; + std::vector<Event> interactions; + for (int i = 0; i < records.size(); i++) + { + interactions.push_back(Event(records[i])); + } + Vector3D apex = getPos(interactions, 56); + EXPECT_DOUBLE_EQ(apex.X, 1.76); + EXPECT_DOUBLE_EQ(apex.Y, 2.62); + EXPECT_DOUBLE_EQ(apex.Z, -0.05); + + Vector3D site2 = getPos(interactions, 55); + EXPECT_DOUBLE_EQ(site2.X, 1.96); + EXPECT_DOUBLE_EQ(site2.Y, 2.46); + EXPECT_DOUBLE_EQ(site2.Z, 0.07); +} + +TEST(Vector3DTest, operators) +{ + Vector3D v(1,2,3); + EXPECT_EQ(v*v, 14); + Vector3D expected(0.1,0.2,0.3); + Vector3D leftMulti = 0.1 * v; + EXPECT_DOUBLE_EQ(expected.X, leftMulti.X); + EXPECT_DOUBLE_EQ(expected.Y, leftMulti.Y); + EXPECT_DOUBLE_EQ(expected.Z, leftMulti.Z); + Vector3D rightMulti = v * 0.1; + EXPECT_DOUBLE_EQ(expected.X, rightMulti.X); + EXPECT_DOUBLE_EQ(expected.Y, rightMulti.Y); + EXPECT_DOUBLE_EQ(expected.Z, rightMulti.Z); + Vector3D division = v / 10; + EXPECT_DOUBLE_EQ(expected.X, division.X); + EXPECT_DOUBLE_EQ(expected.Y, division.Y); + EXPECT_DOUBLE_EQ(expected.Z, division.Z); +} \ No newline at end of file diff --git a/Postprocessing/format.txt b/Postprocessing/format.txt new file mode 100644 index 0000000000000000000000000000000000000000..b6e3f01137a9cdda9f2fd662de7fdd9607b01e57 --- /dev/null +++ b/Postprocessing/format.txt @@ -0,0 +1,103 @@ +dum1 file: + +history_No pid gamma reac target cell E_deposited time x y z weight / nsca / erg_par_in(MeV) + 21 1 2 1 55 35 0.228927 0.04 0.84 2.26 -0.16 1.000E+00 0 0 0 6.610E-01 + 21 1 2 1 55 35 0.199572 0.04 1.11 2.57 -0.04 1.000E+00 0 1 0 4.220E-01 + 21 1 2 3 55 35 0.186618 0.04 1.42 2.27 -0.11 1.000E+00 0 2 0 2.224E-01 + 21 1 2 3 53 35 0.035813 0.04 1.42 2.27 -0.11 1.000E+00 0 2 0 3.581E-02 + 21 2 2 3 53 35 0.010071 0.04 0.84 2.26 -0.16 1.000E+00 0 1 0 1.007E-02 + 48 1 2 1 55 56 0.429851 0.04 1.76 2.62 -0.05 1.000E+00 0 0 0 6.610E-01 + 48 1 2 1 55 55 0.073002 0.04 1.95 2.45 0.11 1.000E+00 0 1 0 2.096E-01 + 48 1 2 3 55 55 0.105975 0.04 1.96 2.46 0.07 1.000E+00 0 2 0 1.366E-01 + 48 1 2 3 55 55 0.030622 0.04 1.97 2.46 0.08 1.000E+00 0 2 0 3.062E-02 + 48 3 2 3 55 56 0.013385 0.04 1.76 2.62 -0.05 1.000E+00 0 1 0 1.338E-02 + 48 2 2 3 55 56 0.008165 0.04 1.76 2.62 -0.05 1.000E+00 0 1 0 8.165E-03 + 51 1 2 1 53 74 0.050400 0.04 2.86 1.86 0.28 1.000E+00 0 1 0 6.281E-01 + 54 1 2 1 55 307 0.358192 0.04 5.11 3.38 0.52 1.000E+00 0 0 0 6.610E-01 + 54 1 2 1 53 98 0.094290 0.05 3.31 3.82 -0.28 1.000E+00 0 1 0 3.028E-01 + 54 1 2 3 53 77 0.179909 0.05 2.91 3.41 0.01 1.000E+00 0 2 0 2.085E-01 + 54 1 2 3 55 77 0.028609 0.05 2.90 3.40 0.00 1.000E+00 0 2 0 2.861E-02 + 59 1 2 1 53 604 0.108711 0.04 7.30 2.09 1.47 1.000E+00 0 0 0 6.610E-01 + 64 1 2 1 55 611 0.099624 0.05 7.21 5.02 -0.33 1.000E+00 0 1 0 6.588E-01 + 67 1 2 1 55 53 0.186788 0.04 2.16 1.40 -0.55 1.000E+00 0 0 0 6.610E-01 + 67 1 2 1 55 53 0.304438 0.04 2.18 1.67 -0.86 1.000E+00 0 1 0 4.742E-01 + + 99999939 2 2 3 53 14 0.001967 0.04 0.23 1.79 -0.98 1.000E+00 0 1 0 1.967E-03 + 99999970 1 2 1 55 33 0.038245 0.04 0.83 1.69 -0.24 1.000E+00 0 0 0 6.610E-01 + 99999970 1 2 3 53 16 0.591504 0.04 0.49 2.65 0.05 1.000E+00 0 1 0 6.228E-01 + 99999970 1 2 3 55 16 0.028609 0.04 0.50 2.65 0.06 1.000E+00 0 1 0 2.861E-02 + 99999970 2 2 3 53 16 0.002642 0.04 0.49 2.65 0.05 1.000E+00 0 1 0 2.642E-03 + 99999988 1 2 1 55 39 0.472426 0.05 1.01 3.91 -1.14 1.000E+00 0 0 0 6.610E-01 + 99999988 1 2 3 55 37 0.155694 0.05 1.06 3.28 -0.93 1.000E+00 0 1 0 1.867E-01 + 99999988 1 2 3 53 37 0.030979 0.05 1.08 3.26 -0.92 1.000E+00 0 1 0 3.098E-02 + 99999988 2 2 3 53 39 0.001901 0.05 1.01 3.91 -1.14 1.000E+00 0 1 0 1.901E-03 + 99999990 1 2 1 55 604 0.202758 0.04 8.00 1.72 -0.13 1.000E+00 0 0 0 6.610E-01 + + 119449 1 2 1 55 78 0.165062 0.05 2.40 3.58 0.68 1.000E+00 0 0 0 6.610E-01 + 119449 2 2 3 53 58 0.034704 0.05 2.40 3.58 0.68 1.000E+00 0 1 0 3.470E-02 + + 26. 0.50000000 2. 119449. 0.030964 3. 2. -999. 0.10000000E+01 58 0.03470000 + 36. 0.50000000 2. 119449. 0.157417 1. 1. -999. 0.10000000E+01 78 0.66100001 + +All_pulses.o + time(ns) pulse_type history_No pulse_height(MeV) reac. weight cell erg_par_in(MeV) + 13. 0.40000000 2. 21. 0.737400 1. 1. -999. 0.10000000E+01 35 0.66100001 + 23. 0.40000000 2. 48. 0.215864 3. 1. -999. 0.10000000E+01 55 0.03062000 + 24. 0.40000000 2. 48. 0.452731 1. 1. -999. 0.10000000E+01 56 0.66100001 + 32. 0.40000000 2. 51. 0.046661 1. 1. -999. 0.10000000E+01 74 0.62809998 + 65. 0.40000000 2. 54. 0.362661 1. 1. -999. 0.10000000E+01 307 0.66100001 + 46. 0.50000000 2. 54. 0.083879 1. 1. -999. 0.10000000E+01 98 0.30280000 + 35. 0.50000000 2. 54. 0.207987 3. 1. -999. 0.10000000E+01 77 0.20850000 + 92. 0.40000000 2. 59. 0.103136 1. 1. -999. 0.10000000E+01 604 0.66100001 + 99. 0.50000000 2. 64. 0.102308 1. 1. -999. 0.10000000E+01 611 0.65880001 + 21. 0.40000000 2. 67. 0.503015 1. 1. -999. 0.10000000E+01 53 0.66100001 + 71. 0.40000000 2. 69. 0.253910 1. 1. -999. 0.10000000E+01 403 0.66100001 + 82. 0.40000000 2. 69. 0.390592 1. 1. -999. 0.10000000E+01 504 0.39390001 + 72. 0.40000000 2. 86. 0.172977 1. 1. -999. 0.10000000E+01 404 0.50120002 + 62. 0.50000000 2. 86. 0.341322 3. 1. -999. 0.10000000E+01 304 0.33640000 + 57. 0.50000000 2. 95. 0.022789 1. 1. -999. 0.10000000E+01 209 0.66100001 + 58. 0.50000000 2. 95. 0.631979 3. 1. -999. 0.10000000E+01 210 0.02831000 + 17. 0.50000000 2. 97. 0.004313 1. 1. -999. 0.10000000E+01 39 0.66100001 + 19. 0.50000000 2. 97. 0.351328 1. 1. -999. 0.10000000E+01 41 0.65170002 + 20. 0.50000000 2. 97. 0.295631 1. 1. -999. 0.10000000E+01 42 0.47999999 + 37. 0.50000000 2. 100. 0.207598 3. 2. -999. 0.10000000E+01 79 0.00400000 + + 31. 0.40000000 2. 89765558. 0.134197 1. 1. -999. 0.10000000E+01 73 0.27219999 + 41. 0.40000000 2. 89765558. 0.385267 1. 1. -999. 0.10000000E+01 93 0.66100001 + 42. 0.40000000 2. 89765558. 0.143121 3. 1. -999. 0.10000000E+01 94 0.13689999 + 79. 0.50000000 2. 89765559. 0.133292 1. 1. -999. 0.10000000E+01 411 0.13230000 + 87. 0.50000000 2. 89765559. 0.540202 1. 1. -999. 0.10000000E+01 509 0.66100001 + 97. 0.50000000 2. 89765569. 0.414940 1. 1. -999. 0.10000000E+01 609 0.66100001 + 21. 0.40000000 2. 89765585. 0.466388 1. 1. -999. 0.10000000E+01 53 0.60430002 + 22. 0.40000000 2. 89765585. 0.114611 3. 1. -999. 0.10000000E+01 54 0.11270000 + 50. 0.50000000 2. 89765590. 0.098971 1. 1. -999. 0.10000000E+01 102 0.66100001 + 92. 0.40000000 2. 89765608. 0.675097 1. 1. -999. 0.10000000E+01 604 0.66100001 + +to reconstruct an event, we need first and second interaction sites, and scattering angle +1. assume initial energy is known as 0.661MeV. we only need to know the energy deposited in the first cell; thus we can select 1-->3 or 1-->1 event pairs. the scattering angle is: + + cos(\theta) = 1 - E_d1 * mc^2 / (0.661 * (0.661-E_d1)) + +2. if the initial energy is unknown, we select 1-->3 event pairs, i.e., photon deposits all its energt in second cell. the scattering angle is: + + cos(\theta) = 1 - E_d1 * mc^2 / ((E_d1+E_d2) * E_d2) + + +while histNo < max: + get pulses with same history number + if number of pulse != 2: + go to next hist number + else: + parse line + sort pulses by erg_par_in + if interaction of second pulse is not 3: + go to next hist number + calculate scattering angle + get interaction sites from dum1 file, based on history number and cell number: + select all interactions with the same history number + sort by cell number + pick up the one with the mathching cell number and max energy of all + +how do we know whether the photon deposits all its energy in the second interaction? +Simulation shows that there are many multiple scattering events. + \ No newline at end of file