Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Standalone spectra out #64

Merged
merged 2 commits into from
Sep 10, 2024
Merged

Standalone spectra out #64

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
9bebd6e
allow separate output for spectra file
KedoKudo Sep 10, 2024
6c73b1e
update unit test
KedoKudo Sep 10, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 5 additions & 0 deletions sophiread/SophireadCLI/include/sophiread_core.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -46,4 +46,9 @@ void updateTOFImages(
std::vector<std::vector<std::vector<unsigned int>>>& tof_images,
const TPX3& batch, double super_resolution,
const std::vector<double>& tof_bin_edges, const std::string& mode);
std::vector<uint64_t> calculateSpectralCounts(
const std::vector<std::vector<std::vector<unsigned int>>>& tof_images);
void writeSpectralFile(const std::string& filename,
const std::vector<uint64_t>& spectral_counts,
const std::vector<double>& tof_bin_edges);
} // namespace sophiread
39 changes: 29 additions & 10 deletions sophiread/SophireadCLI/src/sophiread.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -37,6 +37,7 @@ struct ProgramOptions {
std::string output_tof_imaging;
std::string tof_filename_base = "tof_image";
std::string tof_mode = "neutron";
std::string spectra_filen = "Spectra";
size_t chunk_size = 5ULL * 1024 * 1024 * 1024; // Default 5GB
bool debug_logging = false;
bool verbose = false;
Expand Down Expand Up @@ -68,6 +69,7 @@ void print_usage(const char* program_name) {
spdlog::info(
" -m <tof_mode> TOF mode: 'hit' or 'neutron' (default: "
"neutron)");
spdlog::info(" -s <spectra_filename> Output filename for spectra");
spdlog::info(" -c <chunk_size> Chunk size in MB (default: 5120)");
spdlog::info(" -d Enable debug logging");
spdlog::info(" -v Enable verbose logging");
Expand All @@ -83,7 +85,7 @@ ProgramOptions parse_arguments(int argc, char* argv[]) {
ProgramOptions options;
int opt;

while ((opt = getopt(argc, argv, "i:H:E:u:T:f:m:c:dv")) != -1) {
while ((opt = getopt(argc, argv, "i:H:E:u:T:f:m:s:c:dv")) != -1) {
switch (opt) {
case 'i':
options.input_tpx3 = optarg;
Expand All @@ -106,6 +108,9 @@ ProgramOptions parse_arguments(int argc, char* argv[]) {
case 'm':
options.tof_mode = optarg;
break;
case 's':
options.spectra_filen = optarg;
break;
case 'c':
options.chunk_size = static_cast<size_t>(std::stoull(optarg)) * 1024 *
1024; // Convert MB to bytes
Expand Down Expand Up @@ -210,7 +215,9 @@ int main(int argc, char* argv[]) {

// Initialize TOF images if needed
std::vector<std::vector<std::vector<unsigned int>>> tof_images;
if (!options.output_tof_imaging.empty()) {
const bool needs_tof_images =
!options.output_tof_imaging.empty() || !options.spectra_filen.empty();
if (needs_tof_images) {
tof_images = sophiread::initializeTOFImages(config->getSuperResolution(),
config->getTOFBinEdges());
}
Expand Down Expand Up @@ -258,7 +265,7 @@ int main(int argc, char* argv[]) {
}

// Update TOF images
if (!options.output_tof_imaging.empty()) {
if (needs_tof_images) {
sophiread::updateTOFImages(
tof_images, batch, config->getSuperResolution(),
config->getTOFBinEdges(), options.tof_mode);
Expand Down Expand Up @@ -293,13 +300,6 @@ int main(int argc, char* argv[]) {
eventsFile.close();
}

// Save TOF images if needed
if (!options.output_tof_imaging.empty()) {
sophiread::timedSaveTOFImagingToTIFF(options.output_tof_imaging,
tof_images, config->getTOFBinEdges(),
options.tof_filename_base);
}

auto end = std::chrono::high_resolution_clock::now();
auto elapsed =
std::chrono::duration_cast<std::chrono::microseconds>(end - start)
Expand All @@ -309,6 +309,25 @@ int main(int argc, char* argv[]) {
spdlog::info("Total hits: {}", totalHits);
spdlog::info("Total neutrons: {}", totalNeutrons);

// Save TOF images if needed
if (!options.output_tof_imaging.empty()) {
spdlog::info("Saving TOF imaging to TIFF: {}",
options.output_tof_imaging);
sophiread::timedSaveTOFImagingToTIFF(options.output_tof_imaging,
tof_images, config->getTOFBinEdges(),
options.tof_filename_base);
}

// Save spectra if needed
if (!options.spectra_filen.empty()) {
spdlog::info("Saving spectra to file: {}", options.spectra_filen);
std::string spectral_filename = options.spectra_filen + ".txt";
std::vector<uint64_t> spectral_counts =
sophiread::calculateSpectralCounts(tof_images);
sophiread::writeSpectralFile(spectral_filename, spectral_counts,
config->getTOFBinEdges());
}

} catch (const std::exception& e) {
spdlog::error("Error: {}", e.what());
return 1;
Expand Down
55 changes: 33 additions & 22 deletions sophiread/SophireadCLI/src/sophiread_core.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -365,10 +365,7 @@ void timedSaveTOFImagingToTIFF(
spdlog::info("Created output directory: {}", out_tof_imaging);
}

// 2. Initialize vector for spectral data
std::vector<uint64_t> spectral_counts(tof_images.size(), 0);

// 3. Iterate through each TOF bin and save TIFF files
// 2. Iterate through each TOF bin and save TIFF files
tbb::parallel_for(
tbb::blocked_range<size_t>(0, tof_images.size()),
[&](const tbb::blocked_range<size_t> &range) {
Expand Down Expand Up @@ -439,39 +436,53 @@ void timedSaveTOFImagingToTIFF(
} else {
spdlog::error("Failed to open TIFF file for writing: {}", filename);
}

// Accumulate counts for spectral file
spectral_counts[bin] = std::accumulate(
accumulated_image.cbegin(), accumulated_image.cend(),
static_cast<uint64_t>(0), [](const auto sum, const auto &row) {
return sum + std::accumulate(row.cbegin(), row.cend(), 0ULL);
});
}
});

// 3. Calculate spectral counts
std::vector<uint64_t> spectral_counts = calculateSpectralCounts(tof_images);
// 4. Write spectral file
std::string spectral_filename =
fmt::format("{}/{}_Spectra.txt", out_tof_imaging, tof_filename_base);
// Write spectral data to file
std::ofstream spectral_file(spectral_filename);
writeSpectralFile(spectral_filename, spectral_counts, tof_bin_edges);

auto end = std::chrono::high_resolution_clock::now();
auto duration =
std::chrono::duration_cast<std::chrono::milliseconds>(end - start);
spdlog::info("TIFF and spectra file writing completed in {} ms",
duration.count());
}

std::vector<uint64_t> calculateSpectralCounts(
const std::vector<std::vector<std::vector<unsigned int>>> &tof_images) {
std::vector<uint64_t> spectral_counts(tof_images.size(), 0);

for (size_t bin = 0; bin < tof_images.size(); ++bin) {
spectral_counts[bin] = std::accumulate(
tof_images[bin].cbegin(), tof_images[bin].cend(),
static_cast<uint64_t>(0), [](const auto sum, const auto &row) {
return sum + std::accumulate(row.cbegin(), row.cend(), 0ULL);
});
}

return spectral_counts;
}

void writeSpectralFile(const std::string &filename,
const std::vector<uint64_t> &spectral_counts,
const std::vector<double> &tof_bin_edges) {
std::ofstream spectral_file(filename);
if (spectral_file.is_open()) {
spectral_file << "shutter_time,counts\n";
for (size_t bin = 0; bin < tof_bin_edges.size() - 1; ++bin) {
spectral_file << tof_bin_edges[bin + 1] << "," << spectral_counts[bin]
<< "\n";
}
spectral_file.close();
spdlog::info("Wrote spectral file: {}", spectral_filename);
spdlog::info("Wrote spectral file: {}", filename);
} else {
spdlog::error("Failed to open spectra file for writing: {}",
spectral_filename);
spdlog::error("Failed to open spectra file for writing: {}", filename);
}

auto end = std::chrono::high_resolution_clock::now();
auto duration =
std::chrono::duration_cast<std::chrono::milliseconds>(end - start);
spdlog::info("TIFF and spectra file writing completed in {} ms",
duration.count());
}

} // namespace sophiread
57 changes: 57 additions & 0 deletions sophiread/SophireadCLI/tests/test_sophiread_core.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -150,6 +150,23 @@ TEST_F(SophireadCoreTest, TimedSaveTOFImagingToTIFF) {
EXPECT_TRUE(std::filesystem::exists("test_tof/test_bin_0002.tiff"));
EXPECT_TRUE(std::filesystem::exists("test_tof/test_bin_0003.tiff"));
EXPECT_TRUE(std::filesystem::exists("test_tof/test_Spectra.txt"));

// Check spectra file contents
std::ifstream spectra_file("test_tof/test_Spectra.txt");
std::string line;
std::vector<std::string> lines;
while (std::getline(spectra_file, line)) {
lines.push_back(line);
}
// close the file
spectra_file.close();

EXPECT_EQ(lines.size(), 4); // Header + 3 data lines
EXPECT_EQ(lines[0], "shutter_time,counts");
EXPECT_EQ(lines[1], "0.1,100"); // 10 * 10 * 1 = 100 counts for each bin
EXPECT_EQ(lines[2], "0.2,100");
EXPECT_EQ(lines[3], "0.3,100");

std::filesystem::remove_all("test_tof");
}

Expand Down Expand Up @@ -180,6 +197,46 @@ TEST_F(SophireadCoreTest, UpdateTOFImages) {
EXPECT_TRUE(updated);
}

TEST_F(SophireadCoreTest, CalculateSpectralCounts) {
std::vector<std::vector<std::vector<unsigned int>>> tof_images(
3, std::vector<std::vector<unsigned int>>(
5, std::vector<unsigned int>(
5, 2))); // 3 bins, 5x5 images, all values 2

auto spectral_counts = sophiread::calculateSpectralCounts(tof_images);

EXPECT_EQ(spectral_counts.size(), 3); // 3 bins
for (const auto& count : spectral_counts) {
EXPECT_EQ(count, 50); // 5 * 5 * 2 = 50 for each bin
}
}

TEST_F(SophireadCoreTest, WriteSpectralFile) {
std::vector<uint64_t> spectral_counts = {10, 20, 30};
std::vector<double> tof_bin_edges = {0.0, 0.1, 0.2, 0.3};
std::string filename = "test_spectra.txt";

sophiread::writeSpectralFile(filename, spectral_counts, tof_bin_edges);

EXPECT_TRUE(std::filesystem::exists(filename));

// Read and check file contents
std::ifstream file(filename);
std::string line;
std::vector<std::string> lines;
while (std::getline(file, line)) {
lines.push_back(line);
}

EXPECT_EQ(lines.size(), 4); // Header + 3 data lines
EXPECT_EQ(lines[0], "shutter_time,counts");
EXPECT_EQ(lines[1], "0.1,10");
EXPECT_EQ(lines[2], "0.2,20");
EXPECT_EQ(lines[3], "0.3,30");

std::filesystem::remove(filename);
}

int main(int argc, char** argv) {
::testing::InitGoogleTest(&argc, argv);
return RUN_ALL_TESTS();
Expand Down
Loading