mx/Jungfraujoch
0
0
Fork 0
Code Issues Pull Requests Actions 3 Packages Projects Releases 37 Wiki Activity

JFJochReceiver: Save information on crystal lattice and spot indexing status to HDF5 file

Browse Source
This commit is contained in:
leonarski_f
2023-05-18 10:46:46 +02:00
parent e3998a7488
commit bb4cee9b2d
12 changed files with 68 additions and 37 deletions
Download Patch File Download Diff File Expand all files Collapse all files

3
common/DiffractionSpot.cpp
View File

@@ -55,7 +55,8 @@ DiffractionSpot::operator SpotToSave() const {
return {
.x = x / static_cast<float>(photons),
.y = y / static_cast<float>(photons),
.intensity = static_cast<float>(photons)
.intensity = static_cast<float>(photons),
.indexed = false
};
}

1
common/SpotToSave.h
View File

@@ -8,6 +8,7 @@ struct SpotToSave {
float x;
float y;
float intensity;
bool indexed;
};
#endif //JUNGFRAUJOCH_SPOTTOSAVE_H

1
frame_serialize/ImageMessage.h
View File

@@ -27,6 +27,7 @@ struct DataMessage {
std::vector<SpotToSave> spots;
std::vector<float> rad_int_profile;
uint64_t indexing_result; // 0 - not tried, 1 - tried and failed, 2 - tried and success
std::vector<float> indexing_lattice;
uint64_t bunch_id;
uint32_t jf_info;

11
frame_serialize/JFJochFrameDeserializer.cpp
View File

@@ -271,17 +271,18 @@ void GetCBORMultidimTypedArray(CBORImage &v, CborValue &value) {
void JFJochFrameDeserializer::GetCBORSpots(CborValue &value) {
size_t array_len = GetCBORArrayLen(value);
if (array_len % 3 != 0)
if (array_len % 4 != 0)
throw JFJochException(JFJochExceptionCategory::CBORError,
"Spot array has to have elements multiple of 3");
"Spot array has to have elements multiple of 4");
CborValue array_value;
cborErr(cbor_value_enter_container(&value, &array_value));
for (int i = 0; i < array_len / 3; i++) {
for (int i = 0; i < array_len / 4; i++) {
SpotToSave s{
.x = GetCBORFloat(array_value),
.y = GetCBORFloat(array_value),
.intensity = GetCBORFloat(array_value)
.intensity = GetCBORFloat(array_value),
.indexed = GetCBORBool(array_value)
};
data_message.spots.push_back(s);
}
@@ -362,6 +363,8 @@ void JFJochFrameDeserializer::ProcessImageMessageUserDataElement(CborValue &valu
GetCBORFloatArray(map_value, data_message.rad_int_profile);
else if (key == "indexing_result")
data_message.indexing_result = GetCBORUInt(map_value);
else if (key == "indexing_lattice")
GetCBORFloatArray(map_value, data_message.indexing_lattice);
else if (key == "jf_info")
data_message.jf_info = GetCBORUInt(map_value) & UINT32_MAX;
else if (key == "storage_cell")

6
frame_serialize/JFJochFrameSerializer.cpp
View File

@@ -179,12 +179,13 @@ inline void CBOR_ENC(CborEncoder &encoder, const char* key, const std::vector<Sp
CborEncoder arrayEncoder;
cborErr(cbor_encode_text_stringz(&encoder, key));
cborErr(cbor_encoder_create_array(&encoder, &arrayEncoder, spots.size() * 3));
cborErr(cbor_encoder_create_array(&encoder, &arrayEncoder, spots.size() * 4));
for (auto spot : spots) {
cborErr(cbor_encode_float(&arrayEncoder, spot.x));
cborErr(cbor_encode_float(&arrayEncoder, spot.y));
cborErr(cbor_encode_float(&arrayEncoder, spot.intensity));
cborErr(cbor_encode_boolean(&arrayEncoder, spot.indexed));
}
cborErr(cbor_encoder_close_container(&encoder, &arrayEncoder));
}
@@ -428,11 +429,12 @@ void JFJochFrameSerializer::SerializeImage(const DataMessage& message) {
message.timestamp_base);
cborErr(cbor_encode_text_stringz(&mapEncoder, "user_data"));
cborErr(cbor_encoder_create_map(&mapEncoder, &userDataMapEncoder, 6));
cborErr(cbor_encoder_create_map(&mapEncoder, &userDataMapEncoder, 7));
CBOR_ENC(userDataMapEncoder, "spots", message.spots);
CBOR_ENC(userDataMapEncoder, "rad_int_profile", message.rad_int_profile);
CBOR_ENC(userDataMapEncoder, "indexing_result", message.indexing_result);
CBOR_ENC(userDataMapEncoder, "indexing_lattice", message.indexing_lattice);
CBOR_ENC(userDataMapEncoder, "bunch_id", message.bunch_id);
CBOR_ENC(userDataMapEncoder, "jf_info", (uint64_t) message.jf_info);
CBOR_ENC(userDataMapEncoder, "storage_cell", (uint64_t) message.storage_cell);

9
image_analysis/CrystalLattice.cpp
View File

@@ -60,4 +60,13 @@ CrystalLattice::CrystalLattice(const Eigen::Matrix3f &m) {
vec[i].y = m(i, 1);
vec[i].z = m(i, 2);
}
}
void CrystalLattice::Save(std::vector<float> &output) {
output.resize(9);
for (int i = 0; i < 3; i++) {
output[3 * i + 0] = vec[i].x;
output[3 * i + 1] = vec[i].y;
output[3 * i + 2] = vec[i].z;
}
}

1
image_analysis/CrystalLattice.h
View File

@@ -21,6 +21,7 @@ public:
CrystalLattice ReciprocalLattice() const;
UnitCell GetUnitCell();
Eigen::Matrix3f GetEigenMatrix() const;
void Save(std::vector<float> &output);
};

40
receiver/JFJochReceiver.cpp
View File

@@ -350,11 +350,11 @@ int64_t JFJochReceiver::FrameTransformationThread() {
message.number = image_number;
message.timestamp_base = 10*1000*1000;
message.exptime_base = 10*1000*1000;
message.indexing_result = 0;
bool send_preview = false;
bool send_bkg_estimate = false;
bool calculate_spots = false;
bool index = false;
if ((preview_publisher != nullptr) && (preview_stride > 0) && (image_number % preview_stride == 0))
send_preview = true;
@@ -365,8 +365,6 @@ int64_t JFJochReceiver::FrameTransformationThread() {
calculate_spots = true;
if (rad_int_mapping)
send_bkg_estimate = true;
if (indexer)
index = true;
}
if (experiment.GetSummation() >= threaded_summation_threshold) {
@@ -441,29 +439,31 @@ int64_t JFJochReceiver::FrameTransformationThread() {
if (calculate_spots) {
spot_finder->GetSpotFinderResults(experiment, GetDataProcessingSettings(), spots);
for (const auto & spot : spots)
for (const auto &spot: spots)
message.spots.push_back(spot);
spot_count.AddElement(image_number, spots.size());
}
if (index) {
std::vector<Coord> recip;
for (const auto &i: spots)
recip.push_back(i.ReciprocalCoord(experiment));
if (indexer) {
std::vector<Coord> recip;
for (const auto &i: spots)
recip.push_back(i.ReciprocalCoord(experiment));
auto indexer_result = indexer->Run(recip);
auto indexer_result = indexer->Run(recip);
if (!indexer_result.empty()) {
message.indexing_result = 2;
indexing_solution.AddElement(image_number, 1);
indexing_solution_per_file.Add(image_number % experiment.GetDataFileCount(), 1);
} else {
message.indexing_result = 1;
indexing_solution.AddElement(image_number, 0);
indexing_solution_per_file.Add(image_number % experiment.GetDataFileCount(), 0);
if (!indexer_result.empty()) {
message.indexing_result = 2;
indexing_solution.AddElement(image_number, 1);
indexing_solution_per_file.Add(image_number % experiment.GetDataFileCount(), 1);
for (int i = 0; i < recip.size(); i++)
message.spots[i].indexed = indexer_result[0].indexed_spots[i];
indexer_result[0].l.Save(message.indexing_lattice);
} else {
message.indexing_result = 1;
indexing_solution.AddElement(image_number, 0);
indexing_solution_per_file.Add(image_number % experiment.GetDataFileCount(), 0);
}
}
} else
message.indexing_result = 0;
}
if (send_bkg_estimate) {
uint16_t rad_int_min_bin = std::floor(

7
tests/CBORTest.cpp
View File

@@ -437,8 +437,8 @@ TEST_CASE("CBORSerialize_Image_Spots", "[CBOR]") {
JFJochFrameSerializer serializer(buffer.data(), buffer.size());
std::vector<SpotToSave> spots;
spots.push_back(SpotToSave{.x = 7, .y = 8, .intensity = 34});
spots.push_back(SpotToSave{.x = 37, .y = 48, .intensity = 123});
spots.push_back(SpotToSave{.x = 7, .y = 8, .intensity = 34, .indexed = false});
spots.push_back(SpotToSave{.x = 37, .y = 48, .intensity = 123, .indexed = true});
std::vector<uint8_t> test(1024);
for (int i = 0; i < test.size(); i++)
@@ -474,10 +474,13 @@ JFJochFrameSerializer serializer(buffer.data(), buffer.size());
REQUIRE(image_array.spots.size() == 2);
REQUIRE(image_array.spots[0].intensity == 34);
REQUIRE(!image_array.spots[0].indexed);
REQUIRE(image_array.spots[1].x == 37);
REQUIRE(image_array.spots[1].y == 48);
REQUIRE(image_array.spots[1].intensity == 123);
REQUIRE(image_array.spots[1].indexed);
}
inline bool CmpString(const char *str1, const std::string& str2) {

17
writer/HDF5DataFile.cpp
View File

@@ -31,6 +31,7 @@ HDF5DataFile::~HDF5DataFile() {
std::vector<hsize_t> dims = {spot_count.size(), max_spots};
result_group->SaveVector("nPeaks", spot_count);
result_group->SaveVector("indexingResult", indexing_result);
result_group->SaveVector("indexingLattice", indexing_lattice, {max_image_number + 1, 9});
if (!rad_int_bin_to_q.empty())
rad_int_group->SaveVector("bin_to_q", rad_int_bin_to_q);
@@ -74,9 +75,12 @@ void HDF5DataFile::CreateFile() {
data_set_spot_y = std::make_unique<HDF5DataSet>(*data_file, "/entry/result/peakYPosRaw",
HDF5DataType(0.0f), data_space_spots, dcpl_spots);
data_set_spot_int = std::make_unique<HDF5DataSet>(*data_file, "/entry/result/TotalIntensity",
data_set_spot_int = std::make_unique<HDF5DataSet>(*data_file, "/entry/result/peakTotalIntensity",
HDF5DataType(0.0f), data_space_spots, dcpl_spots);
data_set_spot_indexed = std::make_unique<HDF5DataSet>(*data_file, "/entry/result/peakIndexed",
HDF5DataType(1, false), data_space_spots, dcpl_spots);
spot_count.resize(1);
indexing_result.resize(1);
bunch_id.resize(1);
@@ -84,6 +88,7 @@ void HDF5DataFile::CreateFile() {
timestamp.resize(1);
storage_cell.resize(1);
exptime.resize(1);
indexing_lattice.resize(9);
}
void HDF5DataFile::Write(const DataMessage &msg, uint64_t image_number) {
@@ -97,6 +102,7 @@ void HDF5DataFile::Write(const DataMessage &msg, uint64_t image_number) {
data_set_spot_x->SetExtent({max_image_number+1, max_spots});
data_set_spot_y->SetExtent({max_image_number+1, max_spots});
data_set_spot_int->SetExtent({max_image_number+1, max_spots});
data_set_spot_indexed->SetExtent({max_image_number+1, max_spots});
spot_count.resize(max_image_number + 1);
indexing_result.resize(max_image_number + 1);
@@ -105,6 +111,7 @@ void HDF5DataFile::Write(const DataMessage &msg, uint64_t image_number) {
timestamp.resize(max_image_number + 1);
exptime.resize(max_image_number + 1);
storage_cell.resize(max_image_number + 1);
indexing_lattice.resize((max_image_number + 1) * 9);
}
nimages++;
@@ -115,14 +122,17 @@ void HDF5DataFile::Write(const DataMessage &msg, uint64_t image_number) {
spot_count[image_number] = cnt;
std::vector<float> spot_x(max_spots), spot_y(max_spots), spot_intensity(max_spots);
std::vector<uint8_t> spot_indexed(max_spots);
for (int i = 0; i < cnt; i++) {
spot_x[i] = msg.spots[i].x;
spot_y[i] = msg.spots[i].y;
spot_intensity[i] = msg.spots[i].intensity;
spot_indexed[i] = msg.spots[i].indexed;
}
data_set_spot_x->WriteVec(spot_x, {image_number, 0}, {1, max_spots});
data_set_spot_y->WriteVec(spot_y, {image_number, 0}, {1, max_spots});
data_set_spot_int->WriteVec(spot_intensity, {image_number, 0}, {1, max_spots});
data_set_spot_indexed->WriteVec(spot_indexed, {image_number, 0}, {1, max_spots});
indexing_result[image_number] = msg.indexing_result;
bunch_id[image_number] = msg.bunch_id;
@@ -131,6 +141,11 @@ void HDF5DataFile::Write(const DataMessage &msg, uint64_t image_number) {
storage_cell[image_number] = msg.storage_cell;
exptime[image_number] = msg.exptime;
if (msg.indexing_lattice.size() == 9) {
for (int i = 0; i < 9; i++)
indexing_lattice[image_number * 9 + i] = msg.indexing_lattice[i];
}
if (!msg.rad_int_profile.empty() && (msg.rad_int_profile.size() == rad_int_bin_to_q.size()))
rad_int_group->SaveVector("img" + std::to_string(image_number), msg.rad_int_profile);
}

2
writer/HDF5DataFile.h
View File

@@ -26,6 +26,7 @@ class HDF5DataFile {
std::unique_ptr<HDF5DataSet> data_set_spot_x = nullptr;
std::unique_ptr<HDF5DataSet> data_set_spot_y = nullptr;
std::unique_ptr<HDF5DataSet> data_set_spot_int = nullptr;
std::unique_ptr<HDF5DataSet> data_set_spot_indexed = nullptr;
std::unique_ptr<HDF5Group> result_group = nullptr;
std::unique_ptr<HDF5Group> rad_int_group = nullptr;
@@ -47,6 +48,7 @@ class HDF5DataFile {
std::vector<float> rad_int_bin_to_q;
std::vector<float> rad_int_file_avg;
std::vector<float> indexing_lattice;
std::vector<uint32_t> spot_count;
const size_t max_spots;

7
writer/HDF5Objects.h
View File

@@ -115,8 +115,6 @@ public:
std::vector<hsize_t> dim = {}, CompressionAlgorithm algorithm = CompressionAlgorithm::NO_COMPRESSION);
};
#include <iostream>
class HDF5Group : public HDF5Object {
public:
HDF5Group(const HDF5Object& parent, const std::string& name);
@@ -157,10 +155,6 @@ public:
const std::vector<hsize_t> &size) {
HDF5DataSpace mem_space({v.size()});
HDF5DataSpace file_space(*this);
std::cout << int32_t(file_space.GetNumOfDimensions()) << std::endl;
std::cout << int32_t(file_space.GetDimensions()[0]) << std::endl;
std::cout << int32_t(mem_space.GetDimensions()[0]) << std::endl;
file_space.SelectHyperslab(start, size);
if (H5Dwrite(id,
@@ -178,7 +172,6 @@ public:
const std::vector<hsize_t> &start) {
HDF5DataSpace mem_space({1});
HDF5DataSpace file_space(*this);
std::cout << uint32_t(file_space.GetNumOfDimensions()) << std::endl;
file_space.SelectHyperslab(start, {1});
if (H5Dwrite(id,
HDF5DataType(val).GetID(),