Commit 2119faee authored by Purves, Murray's avatar Purves, Murray
Browse files

WIP PCF file reader

parent f7b01e7a
Pipeline #16089 passed with stages
in 8 minutes and 38 seconds
#include "radixio/spectrum.hh"
namespace radix
{
short Spectrum::nrps() const { return mNrps; }
void Spectrum::setNrps(short nrps) { mNrps = nrps; }
bool Spectrum::dhsVersion() const { return mDhsVersion; }
void Spectrum::setDhsVersion(bool dhsType) { mDhsVersion = dhsType; }
std::string Spectrum::energyCalibrationLabel() const
{
return mEnergyCalibrationLabel;
}
void Spectrum::setEnergyCalibrationLabel(
const std::string &energyCalibrationLabel)
{
mEnergyCalibrationLabel = energyCalibrationLabel;
}
float Spectrum::energyCalibrationOffset() const
{
return mEnergyCalibrationOffset;
}
void Spectrum::setEnergyCalibrationOffset(float energyCalibrationOffset)
{
mEnergyCalibrationOffset = energyCalibrationOffset;
}
float Spectrum::energyCalibrationGain() const { return mEnergyCalibrationGain; }
void Spectrum::setEnergyCalibrationGain(float energyCalibrationGain)
{
mEnergyCalibrationGain = energyCalibrationGain;
}
float Spectrum::energyCalibrationQuadraticTerm() const
{
return mEnergyCalibrationQuadraticTerm;
}
void Spectrum::setEnergyCalibrationQuadraticTerm(
float energyCalibrationQuadraticTerm)
{
mEnergyCalibrationQuadraticTerm = energyCalibrationQuadraticTerm;
}
float Spectrum::energyCalibrationCubicTerm() const
{
return mEnergyCalibrationCubicTerm;
}
void Spectrum::setEnergyCalibrationCubicTerm(float energyCalibrationCubicTerm)
{
mEnergyCalibrationCubicTerm = energyCalibrationCubicTerm;
}
float Spectrum::energyCalibrationLowEnergy() const
{
return mEnergyCalibrationLowEnergy;
}
void Spectrum::setEnergyCalibrationLowEnergy(float energyCalibrationLowEnergy)
{
mEnergyCalibrationLowEnergy = energyCalibrationLowEnergy;
}
std::string Spectrum::lastModifiedHash() const { return mLastModifiedHash; }
void Spectrum::setLastModifiedHash(const std::string &lastModifiedHash)
{
mLastModifiedHash = lastModifiedHash;
}
std::string Spectrum::uuid() const { return mUuid; }
void Spectrum::setUuid(const std::string &uuid) { mUuid = uuid; }
std::string Spectrum::inspection() const { return mInspection; }
void Spectrum::setInspection(const std::string &inspection)
{
mInspection = inspection;
}
short Spectrum::laneNumber() const { return mLaneNumber; }
void Spectrum::setLaneNumber(short laneNumber) { mLaneNumber = laneNumber; }
std::string Spectrum::measurementRemark() const { return mMeasurementRemark; }
void Spectrum::setMeasurementRemark(const std::string &measurementRemark)
{
mMeasurementRemark = measurementRemark;
}
std::string Spectrum::instrumentType() const { return mInstrumentType; }
void Spectrum::setInstrumentType(const std::string &instrumentType)
{
mInstrumentType = instrumentType;
}
std::string Spectrum::manufacturer() const { return mManufacturer; }
void Spectrum::setManufacturer(const std::string &manufacturer)
{
mManufacturer = manufacturer;
}
std::string Spectrum::instrumentModel() const { return mInstrumentModel; }
void Spectrum::setInstrumentModel(const std::string &instrumentModel)
{
mInstrumentModel = instrumentModel;
}
std::string Spectrum::instrumentID() const { return mInstrumentID; }
void Spectrum::setInstrumentID(const std::string &instrumentID)
{
mInstrumentID = instrumentID;
}
std::string Spectrum::itemDescription() const { return mItemDescription; }
void Spectrum::setItemDescription(const std::string &itemDescription)
{
mItemDescription = itemDescription;
}
std::string Spectrum::measurementLocationName() const
{
return mMeasurementLocationName;
}
void Spectrum::setMeasurementLocationName(
const std::string &measurementLocationName)
{
mMeasurementLocationName = measurementLocationName;
}
std::string Spectrum::measurementLocationCoords() const
{
return mMeasurementLocationCoords;
}
void Spectrum::setMeasurementLocationCoords(
const std::string &measurementLocationCoords)
{
mMeasurementLocationCoords = measurementLocationCoords;
}
short Spectrum::itemDetectorDistance() const { return mItemDetectorDistance; }
void Spectrum::setItemDetectorDistance(short itemDetectorDistance)
{
mItemDetectorDistance = itemDetectorDistance;
}
short Spectrum::occupancyNumber() const { return mOccupancyNumber; }
void Spectrum::setOccupancyNumber(short occupancyNumber)
{
mOccupancyNumber = occupancyNumber;
}
std::string Spectrum::cargoType() const { return mCargoType; }
void Spectrum::setCargoType(const std::string &cargoType)
{
mCargoType = cargoType;
}
std::string Spectrum::deviationPairPresence() const
{
return mDeviationPairPresence;
}
void Spectrum::setDeviationPairPresence(
const std::string &deviationPairPresence)
{
mDeviationPairPresence = deviationPairPresence;
}
std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
Spectrum::deviationPairEnergies() const
{
return mDeviationPairEnergies;
}
void Spectrum::setDeviationPairEnergies(
const std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
&deviationPairEnergies)
{
mDeviationPairEnergies = deviationPairEnergies;
}
std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
Spectrum::deviationPairOffsets() const
{
return mDeviationPairOffsets;
}
void Spectrum::setDeviationPairOffsets(
const std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
&deviationPairOffsets)
{
mDeviationPairOffsets = deviationPairOffsets;
}
short Spectrum::srsi() const { return mSrsi; }
void Spectrum::setSrsi(short srsi) { mSrsi = srsi; }
int Spectrum::spectrumDataCount() { return mSpectrumData.size(); }
void Spectrum::spectrumData(int i, std::string &buffer, std::string &dateTime,
std::string &tag, float &liveTime, float &totalTime,
float &energyCalOffset, float &energyCalGain,
float &energyCalQuadraticTerm,
float &energyCalCubicTerm,
float &energyCalLowEnergyTerm, float &occupancyFlag,
float &totalNeutronCount, int &numberOfChannels,
std::vector<float> countsByChannel)
{
if (i > 0 && i < mSpectrumData.size())
{
size_t i_ = size_t(i);
SpectrumData spectrum = mSpectrumData[i];
buffer = spectrum.buffer;
dateTime = spectrum.dateTime;
tag = spectrum.tag;
liveTime = spectrum.liveTime;
totalTime = spectrum.totalTime;
energyCalOffset = spectrum.energyCalOffset;
energyCalGain = spectrum.energyCalGain;
energyCalQuadraticTerm = spectrum.energyCalQuadraticTerm;
energyCalCubicTerm = spectrum.energyCalCubicTerm;
energyCalLowEnergyTerm = spectrum.energyCalLowEnergyTerm;
occupancyFlag = spectrum.occupancyFlag;
totalNeutronCount = spectrum.totalNeutronCount;
numberOfChannels = spectrum.numberOfChannels;
countsByChannel = spectrum.countsByChannel;
}
}
void Spectrum::addSpectrumData(
std::string buffer, std::string dateTime, std::string tag, float liveTime,
float totalTime, float energyCalOffset, float energyCalGain,
float energyCalQuadraticTerm, float energyCalCubicTerm,
float energyCalLowEnergyTerm, float occupancyFlag, float totalNeutronCount,
int numberOfChannels, std::vector<float> countsByChannel)
{
SpectrumData spectrumData;
spectrumData.buffer = buffer;
spectrumData.dateTime = dateTime;
spectrumData.tag = tag;
spectrumData.liveTime = liveTime;
spectrumData.totalTime = totalTime;
spectrumData.energyCalOffset = energyCalOffset;
spectrumData.energyCalGain = energyCalGain;
spectrumData.energyCalQuadraticTerm = energyCalQuadraticTerm;
spectrumData.energyCalCubicTerm = energyCalCubicTerm;
spectrumData.energyCalLowEnergyTerm = energyCalLowEnergyTerm;
spectrumData.occupancyFlag = occupancyFlag;
spectrumData.totalNeutronCount = totalNeutronCount;
spectrumData.numberOfChannels = numberOfChannels;
spectrumData.countsByChannel = countsByChannel;
mSpectrumData.push_back(spectrumData);
}
}
#ifndef RADIX_RADIXIO_SPECTRUM_HH_
#define RADIX_RADIXIO_SPECTRUM_HH_
#include <array>
#include <memory>
#include <string>
#include <vector>
#include "radixcore/visibility.hh"
......@@ -13,6 +15,137 @@ class RADIX_PUBLIC Spectrum
public:
typedef std::shared_ptr<Spectrum> SP;
struct SpectrumData
{
std::string buffer, dateTime, tag;
float liveTime, totalTime, energyCalOffset, energyCalGain,
energyCalQuadraticTerm, energyCalCubicTerm, energyCalLowEnergyTerm,
occupancyFlag, totalNeutronCount;
int numberOfChannels;
std::vector<float> countsByChannel, l;
}; // struct SpectrumData
short nrps() const;
void setNrps(short nrps);
bool dhsVersion() const;
void setDhsVersion(bool dhsVersion);
std::string energyCalibrationLabel() const;
void setEnergyCalibrationLabel(const std::string &energyCalibrationLabel);
float energyCalibrationOffset() const;
void setEnergyCalibrationOffset(float energyCalibrationOffset);
float energyCalibrationGain() const;
void setEnergyCalibrationGain(float energyCalibrationGain);
float energyCalibrationQuadraticTerm() const;
void setEnergyCalibrationQuadraticTerm(float energyCalibrationQuadraticTerm);
float energyCalibrationCubicTerm() const;
void setEnergyCalibrationCubicTerm(float energyCalibrationCubicTerm);
float energyCalibrationLowEnergy() const;
void setEnergyCalibrationLowEnergy(float energyCalibrationLowEnergy);
std::string lastModifiedHash() const;
void setLastModifiedHash(const std::string &lastModifiedHash);
std::string uuid() const;
void setUuid(const std::string &uuid);
std::string inspection() const;
void setInspection(const std::string &inspection);
short laneNumber() const;
void setLaneNumber(short laneNumber);
std::string measurementRemark() const;
void setMeasurementRemark(const std::string &measurementRemark);
std::string instrumentType() const;
void setInstrumentType(const std::string &instrumentType);
std::string manufacturer() const;
void setManufacturer(const std::string &manufacturer);
std::string instrumentModel() const;
void setInstrumentModel(const std::string &instrumentModel);
std::string instrumentID() const;
void setInstrumentID(const std::string &instrumentID);
std::string itemDescription() const;
void setItemDescription(const std::string &itemDescription);
std::string measurementLocationName() const;
void setMeasurementLocationName(const std::string &measurementLocationName);
std::string measurementLocationCoords() const;
void setMeasurementLocationCoords(
const std::string &measurementLocationCoords);
short itemDetectorDistance() const;
void setItemDetectorDistance(short itemDetectorDistance);
short occupancyNumber() const;
void setOccupancyNumber(short occupancyNumber);
std::string cargoType() const;
void setCargoType(const std::string &cargoType);
std::string deviationPairPresence() const;
void setDeviationPairPresence(const std::string &deviationPairPresence);
std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
deviationPairEnergies() const;
void setDeviationPairEnergies(
const std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
&deviationPairEnergies);
std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
deviationPairOffsets() const;
void setDeviationPairOffsets(
const std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
&deviationPairOffsets);
short srsi() const;
void setSrsi(short srsi);
int spectrumDataCount();
void spectrumData(int i, std::string &buffer, std::string &dateTime,
std::string &tag, float &liveTime, float &totalTime,
float &energyCalOffset, float &energyCalGain,
float &energyCalQuadraticTerm, float &energyCalCubicTerm,
float &energyCalLowEnergyTerm, float &occupancyFlag,
float &totalNeutronCount, int &numberOfChannels,
std::vector<float> countsByChannel);
void addSpectrumData(std::string buffer, std::string dateTime,
std::string tag, float liveTime, float totalTime,
float energyCalOffset, float energyCalGain,
float energyCalQuadraticTerm, float energyCalCubicTerm,
float energyCalLowEnergyTerm, float occupancyFlag,
float totalNeutronCount, int numberOfChannels,
std::vector<float> countsByChannel);
private:
// PCF file objects -
// TODO work out what/how to change to support SPE and PCF
// Number of records per spectrum
short mNrps;
// Is the file a DHS version?
bool mDhsVersion = false;
// Non-DHS header objects
std::string mEnergyCalibrationLabel;
float mEnergyCalibrationOffset;
float mEnergyCalibrationGain;
float mEnergyCalibrationQuadraticTerm;
float mEnergyCalibrationCubicTerm;
float mEnergyCalibrationLowEnergy;
// DHS header objects
std::string mLastModifiedHash;
std::string mUuid;
std::string mInspection;
short mLaneNumber;
std::string mMeasurementRemark;
std::string mInstrumentType;
std::string mManufacturer;
std::string mInstrumentModel;
std::string mInstrumentID;
std::string mItemDescription;
std::string mMeasurementLocationName;
std::string mMeasurementLocationCoords;
short mItemDetectorDistance;
short mOccupancyNumber;
std::string mCargoType;
// Deviation pair presence indicator
std::string mDeviationPairPresence;
// Deviation pairs
std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
mDeviationPairEnergies;
std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
mDeviationPairOffsets;
// Spectral record start index
short mSrsi;
// Spectrum data
std::vector<SpectrumData> mSpectrumData;
}; // class RADIX_PUBLIC Spectrum
template <typename data_type>
......
......@@ -3,7 +3,10 @@
#include <string>
#include "radixbug/bug.hh"
#include "radixcore/stringfunctions.hh"
#include "radixcore/visibility.hh"
#include "radixio/eafstream.hh"
#include "radixio/spectrum.hh"
namespace radix
......@@ -19,6 +22,267 @@ bool SpectrumPCFStream<data_type>::read_from(const std::string &file)
{
bool result = false;
radix_line("Reading data from PCF file: " << file);
// Open a stream of the file
eafstream stream(file.c_str(), std::ios::in | std::ios::binary);
if (stream.is_open() == false)
{
return result;
}
// File should be little endian
stream.setReverseBytes(false);
// Read the header
radix_line(" Reading header...");
// Number of records per spectrum (NRPS)
short nrps = stream.readShort();
mData->setNrps(nrps);
radix_line(" nrps: " << nrps);
// Version (3-character string)
std::string version = stream.readString(3);
radix_line(" version: " << version);
if (version.compare("DHS") == 0)
{
radix_line(" File has DHS-type header");
mData->setDhsVersion(true);
// Read in the DHS type header data
std::string lastModifiedHash = stream.readString(7);
mData->setLastModifiedHash(lastModifiedHash);
std::string uuid = stream.readString(36);
mData->setUuid(uuid);
std::string inspection = stream.readString(16);
mData->setInspection(inspection);
short laneNumber = stream.readShort();
mData->setLaneNumber(laneNumber);
std::string measurementRemark = stream.readString(26);
mData->setMeasurementRemark(measurementRemark);
std::string instrumentType = stream.readString(28);
mData->setInstrumentType(instrumentType);
std::string manufacturer = stream.readString(28);
mData->setManufacturer(manufacturer);
std::string instrumentModel = stream.readString(18);
mData->setInstrumentModel(instrumentModel);
std::string instrumentID = stream.readString(18);
mData->setInstrumentID(instrumentID);
std::string itemDescription = stream.readString(20);
mData->setItemDescription(itemDescription);
std::string measurementLocationName = stream.readString(16);
mData->setMeasurementLocationName(measurementLocationName);
std::string measurementLocationCoords = stream.readString(16);
mData->setMeasurementLocationCoords(measurementLocationCoords);
short itemDetectorDistance = stream.readShort();
mData->setItemDetectorDistance(itemDetectorDistance);
short occupancyNumber = stream.readShort();
mData->setOccupancyNumber(occupancyNumber);
std::string cargoType = stream.readString(16);
mData->setCargoType(cargoType);
}
else
{
radix_line(" File has non-DHS-type header");
mData->setDhsVersion(false);
// Read in the non-DHS type header data
std::string mEnergyCalibrationLabel = stream.readString(4);
mData->setEnergyCalibrationLabel(mEnergyCalibrationLabel);
float mEnergyCalibrationOffset = stream.readFloat();
mData->setEnergyCalibrationOffset(mEnergyCalibrationOffset);
float mEnergyCalibrationGain = stream.readFloat();
mData->setEnergyCalibrationGain(mEnergyCalibrationGain);
float mEnergyCalibrationQuadraticTerm = stream.readFloat();
mData->setEnergyCalibrationQuadraticTerm(mEnergyCalibrationQuadraticTerm);
float mEnergyCalibrationCubicTerm = stream.readFloat();
mData->setEnergyCalibrationCubicTerm(mEnergyCalibrationCubicTerm);
float mEnergyCalibrationLowEnergy = stream.readFloat();
mData->setEnergyCalibrationLowEnergy(mEnergyCalibrationLowEnergy);
// Skip to end of header section (at 256 bytes)
radix_line(" Skipping " << 256 - stream.bytesRead()
<< " bytes to end of header (@256 bytes)");
stream.skipBytes(256 - stream.bytesRead());
}
radix_line(" Header read complete; looking for deviation pairs...");
// Read the deviation pairs (if present)
std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
deviationPairEnergies = {0};
std::array<std::array<std::array<std::array<float, 4>, 8>, 8>, 20>
deviationPairOffsets = {0};
std::string deviationPairPresence = radix::trim_string(stream.readString(30));
if (deviationPairPresence.compare("DeviationPairsInFile") == 0)
{
radix_line(" Deviation pairs present; reading...");
// Set SRSI to correct value
mData->setSrsi(83);
// Read uncompressed deviation pairs
// Move to beginning of pairs (at 512 bytes)
radix_line(" Skipping " << 512 - stream.bytesRead()
<< " bytes to beginning of pairs (@512 bytes)");
stream.skipBytes(512 - stream.bytesRead());
for (size_t column = 0; column < 2; ++column)
{
radix(" ");
for (size_t panel = 0; panel < 8; ++panel)
{
for (size_t mca = 0; mca < 8; ++mca)
{
for (size_t pair = 0; pair < 20; ++pair)
{
float energy = stream.readFloat();
float offset = stream.readFloat();
deviationPairEnergies[column][panel][mca][pair] = energy;
deviationPairOffsets[column][panel][mca][pair] = offset;
}
radix(".");
}
radix("-");
}
radix_line("|");
}
}
else if (deviationPairPresence.compare("DeviationPairsInFileCompressed") == 0)
{
radix_line(" Compressed deviation pairs present; reading...");
// Set SRSI to correct value
mData->setSrsi(83);
// Read compressed deviation pairs
// Move to beginning of pairs
stream.skipBytes(226);
for (size_t column = 0; column < 4; ++column)
{
radix(" ");
for (size_t panel = 0; panel < 8; ++panel)
{
for (size_t mca = 0; mca < 8; ++mca)
{
for (size_t pair = 0; pair < 20; ++pair)
{
float energy = float(stream.readShort()) / 10.f;
short offset = float(stream.readShort()) / 10.f;
deviationPairEnergies[column][panel][mca][pair] = energy;
deviationPairOffsets[column][panel][mca][pair] = offset;
}
radix(".");
}
radix("-");
}
radix_line("|");
}
}
else
{
radix_line(" No deviation pairs present");
// Set SRSI to correct value
mData->setSrsi(2);
}
// Set deviation pairs
mData->setDeviationPairEnergies(deviationPairEnergies);
mData->setDeviationPairOffsets(deviationPairOffsets);
radix_line(" Read deviation pairs; reading spectral data...");
// Read spectral header
// Make sure we only read till EOF
stream.peek();
while (stream.good())
{
// Read spectral header
radix_line(" Reading spectral header...");
// Buffer