Line data    Source code

       1             : #include "EDFHandler.h"
       2             : #include "Logger.h"
       3             : #include "VolumeDescriptor.h"
       4             : 
       5             : #include <stdexcept>
       6             : 
       7             : namespace elsa
       8             : {
       9             :     template <typename data_t>
      10             :     DataContainer<data_t> EDF::read(std::string filename)
      11           4 :     {
      12           4 :         Logger::get("EDF")->info("Reading data from {}", filename);
      13             : 
      14             :         // open the file
      15           4 :         std::ifstream file(filename, std::ios::binary);
      16           4 :         if (!file.good())
      17           4 :             throw Error("EDF::read: cannot read from '" + filename + "'");
      18             : 
      19           0 :         return EDF::read<data_t>(file);
      20           0 :     }
      21             : 
      22             :     template <typename data_t>
      23             :     DataContainer<data_t> EDF::read(std::istream& input)
      24           3 :     {
      25             :         // get the meta data from the header
      26           3 :         auto properties = readHeader(input);
      27           3 :         auto [descriptor, dataType] = parseHeader(properties);
      28             : 
      29             :         // read in the data
      30           3 :         DataContainer<data_t> dataContainer(*descriptor);
      31             : 
      32           3 :         if (dataType == DataUtils::DataType::UINT16)
      33           0 :             DataUtils::parseRawData<uint16_t, data_t>(input, dataContainer);
      34           3 :         else if (dataType == DataUtils::DataType::FLOAT32)
      35           3 :             DataUtils::parseRawData<float, data_t>(input, dataContainer);
      36           0 :         else if (dataType == DataUtils::DataType::FLOAT64)
      37           0 :             DataUtils::parseRawData<double, data_t>(input, dataContainer);
      38           0 :         else
      39           0 :             throw Error("EDF::read: invalid/unsupported data type");
      40             : 
      41           3 :         return dataContainer;
      42           3 :     }
      43             : 
      44             :     template <typename data_t>
      45             :     void EDF::write(const DataContainer<data_t>& data, std::string filename)
      46           2 :     {
      47           2 :         Logger::get("EDF")->info("Writing data to {}", filename);
      48             : 
      49             :         // open the file
      50           2 :         std::ofstream file(filename, std::ios::binary);
      51           2 :         if (!file.good())
      52           0 :             throw Error("EDF::write: cannot write to '" + filename + "'");
      53             : 
      54           2 :         EDF::write<data_t>(data, file);
      55           2 :     }
      56             : 
      57             :     template <typename data_t>
      58             :     void EDF::write(const DataContainer<data_t>& data, std::ostream& output)
      59           7 :     {
      60             :         // output the header
      61           7 :         writeHeader(output, data);
      62             : 
      63             :         // output the raw data
      64             :         // TODO: this would be more efficient if we had a data pointer...
      65        2352 :         for (index_t i = 0; i < data.getSize(); ++i)
      66        2345 :             output.write(reinterpret_cast<const char*>(&data[i]), sizeof(data_t));
      67           7 :     }
      68             : 
      69             :     std::map<std::string, std::string> EDF::readHeader(std::istream& file)
      70           3 :     {
      71           3 :         std::map<std::string, std::string> properties;
      72             : 
      73             :         // read a single character and make sure that a header is opened
      74           3 :         if (file.eof() || file.get() != '{')
      75           0 :             throw InvalidArgumentError("EDF::readHeader: no header opening marker");
      76             : 
      77             :         // read header data
      78          26 :         while (!file.eof()) {
      79             :             // skip whitespace
      80        2688 :             while (!file.eof()) {
      81        2688 :                 const int chr = file.peek();
      82        2688 :                 if (chr != '\r' && chr != '\n' && chr != ' ' && chr != '\t')
      83          26 :                     break;
      84        2662 :                 file.ignore(1);
      85        2662 :             }
      86             : 
      87             :             // abort if the header is closed
      88          26 :             if (file.eof() || file.peek() == '}')
      89           3 :                 break;
      90             : 
      91             :             // extract the property assignment
      92          23 :             bool quotesSingle = false, quotesDouble = false;
      93          23 :             std::string assignment;
      94         401 :             while (!file.eof()) {
      95         401 :                 auto chr = static_cast<char>(file.get());
      96             : 
      97             :                 // abort on end-of-assignment
      98         401 :                 if (chr == ';' && !(quotesSingle || quotesDouble))
      99          23 :                     break;
     100             : 
     101             :                 // check for quote characters
     102         378 :                 if (chr == '\'')
     103           0 :                     quotesSingle = !quotesSingle;
     104         378 :                 if (chr == '\"')
     105           0 :                     quotesDouble = !quotesDouble;
     106             : 
     107         378 :                 assignment += chr;
     108         378 :             }
     109             : 
     110             :             // split the assignment
     111          23 :             auto delim = assignment.find('=');
     112          23 :             if (delim == std::string::npos)
     113           0 :                 throw InvalidArgumentError("failed reading name/value delimiter");
     114             : 
     115          23 :             std::string name = assignment.substr(0, delim);
     116          23 :             StringUtils::trim(name);
     117             : 
     118          23 :             std::string value = assignment.substr(delim + 1);
     119          23 :             StringUtils::trim(value);
     120             : 
     121             :             // remove quotes (if they exist)
     122          23 :             if (value[0] == value[value.size() - 1] && (value[0] == '\'' || value[0] == '\"'))
     123           0 :                 value = value.substr(1, value.size() - 2);
     124             : 
     125          23 :             StringUtils::toLower(name);
     126          23 :             properties[name] = value;
     127          23 :         }
     128           3 :         file.ignore(2); // end of header marker
     129             : 
     130           3 :         return properties;
     131           3 :     }
     132             : 
     133             :     std::pair<std::unique_ptr<DataDescriptor>, DataUtils::DataType>
     134             :         EDF::parseHeader(const std::map<std::string, std::string>& properties)
     135           3 :     {
     136             :         // read the dimensions
     137           3 :         std::vector<index_t> dim;
     138           8 :         for (index_t i = 1;; i++) {
     139             :             // assemble the property name
     140           8 :             std::stringstream aux;
     141           8 :             aux << "dim_" << i;
     142             : 
     143             :             // try to find the property
     144           8 :             auto dimIt = properties.find(aux.str());
     145           8 :             if (dimIt == properties.end())
     146           3 :                 break;
     147             : 
     148           5 :             dim.push_back(DataUtils::parse<index_t>(dimIt->second));
     149           5 :         }
     150           3 :         const auto nDims = static_cast<index_t>(dim.size());
     151           3 :         if (nDims == 0u)
     152           0 :             throw Error("EDF::parseHeader: dimension information not found");
     153             : 
     154             :         // parse the (non-standard) spacing tag
     155           3 :         std::vector<real_t> spacing;
     156           3 :         auto spacingIt = properties.find("spacing");
     157           3 :         if (spacingIt != properties.end())
     158           3 :             spacing = DataUtils::parseVector<real_t>(spacingIt->second);
     159             : 
     160             :         // check for a byte order tag, but fall back to the default value
     161           3 :         auto byteorderIt = properties.find("byteorder");
     162           3 :         if (byteorderIt != properties.end()) {
     163           3 :             std::string byteorderValue = byteorderIt->second;
     164           3 :             StringUtils::toLower(byteorderValue);
     165             : 
     166           3 :             if (byteorderValue != "lowbytefirst")
     167           0 :                 throw Error("EDF::parseHeader: unsupported byte order value");
     168           3 :         }
     169             : 
     170             :         // check for the 'element type' value
     171           3 :         DataUtils::DataType dataType;
     172           3 :         auto datatypeIt = properties.find("datatype");
     173           3 :         if (datatypeIt != properties.end()) {
     174           3 :             std::string datatypeValue = datatypeIt->second;
     175           3 :             StringUtils::toLower(datatypeValue);
     176             : 
     177           3 :             if (datatypeValue == "signedbyte")
     178           0 :                 dataType = DataUtils::DataType::INT8;
     179           3 :             else if (datatypeValue == "unsignedbyte")
     180           0 :                 dataType = DataUtils::DataType::UINT8;
     181           3 :             else if (datatypeValue == "signedshort")
     182           0 :                 dataType = DataUtils::DataType::INT16;
     183           3 :             else if (datatypeValue == "unsignedshort")
     184           0 :                 dataType = DataUtils::DataType::UINT16;
     185           3 :             else if (datatypeValue == "float" || datatypeValue == "floatvalue"
     186           3 :                      || datatypeValue == "real")
     187           3 :                 dataType = DataUtils::DataType::FLOAT32;
     188           0 :             else if (datatypeValue == "double" || datatypeValue == "doublevalue")
     189           0 :                 dataType = DataUtils::DataType::FLOAT64;
     190           0 :             else
     191           0 :                 throw Error("EDF::parseHeader: invalid/unsupported data type");
     192           0 :         } else
     193           0 :             throw Error("EDF::parseHeader: data type not found");
     194             : 
     195           3 :         auto compressionIt = properties.find("compression");
     196           3 :         if (compressionIt != properties.end())
     197           0 :             throw Error("EDF::parseHeader: compression not supported");
     198             : 
     199           3 :         index_t size = 0;
     200           3 :         auto sizeIt = properties.find("size");
     201           3 :         if (sizeIt != properties.end())
     202           3 :             size = DataUtils::parse<index_t>(sizeIt->second);
     203             : 
     204           3 :         auto imageIt = properties.find("image");
     205           3 :         if (imageIt != properties.end() && DataUtils::parse<index_t>(imageIt->second) != 1)
     206           0 :             throw Error("EDF::parseHeader: image not set to 1");
     207             : 
     208             :         // convert size
     209           3 :         IndexVector_t dimSizeVec(nDims);
     210           8 :         for (index_t i = 0; i < nDims; ++i)
     211           5 :             dimSizeVec[i] = dim[static_cast<std::size_t>(i)];
     212           3 :         if (dimSizeVec.prod() * DataUtils::getSizeOfDataType(dataType) != size)
     213           0 :             throw Error("EDF::parseHeader: size inconsistency");
     214             : 
     215             :         // convert spacing
     216           3 :         RealVector_t dimSpacingVec(RealVector_t::Ones(nDims));
     217           3 :         if (!spacing.empty()) {
     218           3 :             if (nDims != static_cast<index_t>(spacing.size()))
     219           0 :                 throw Error("EDF::parseHeader: spacing inconsistency");
     220           8 :             for (index_t i = 0; i < nDims; ++i)
     221           5 :                 dimSpacingVec[i] = spacing[static_cast<std::size_t>(i)];
     222           3 :         }
     223             : 
     224           3 :         return std::make_pair(std::make_unique<VolumeDescriptor>(dimSizeVec, dimSpacingVec),
     225           3 :                               dataType);
     226           3 :     }
     227             : 
     228             :     template <typename data_t>
     229             :     void EDF::writeHeader(std::ostream& file, const DataContainer<data_t>& data)
     230           7 :     {
     231             :         // open the header
     232           7 :         file << "{\n";
     233             : 
     234           7 :         file << "HeaderID = EH:000001:000000:000000;\n";
     235           7 :         file << "Image = " << 1 << ";\n";
     236           7 :         file << "ByteOrder = LowByteFirst;\n";
     237           7 :         file << "DataType = " << getDataTypeName(data) << ";\n";
     238             : 
     239           7 :         auto& descriptor = data.getDataDescriptor();
     240             : 
     241             :         // write dimension and size
     242          19 :         for (index_t i = 0; i < descriptor.getNumberOfDimensions(); ++i)
     243          12 :             file << "Dim_" << (i + 1) << " = "
     244          12 :                  << descriptor.getNumberOfCoefficientsPerDimension()[i] << ";\n";
     245           7 :         file << "Size = "
     246           7 :              << descriptor.getNumberOfCoefficients() * static_cast<index_t>(sizeof(data_t))
     247           7 :              << ";\n";
     248             : 
     249             :         // write spacing
     250           7 :         file << "Spacing =";
     251          19 :         for (index_t i = 0; i < descriptor.getNumberOfDimensions(); ++i)
     252          12 :             file << ' ' << descriptor.getSpacingPerDimension()[i];
     253           7 :         file << ";\n";
     254             : 
     255             :         // pad the header by adding spaces such that the header ends on a kilobyte boundary
     256           7 :         index_t n = 1024;
     257           7 :         while (n < (static_cast<index_t>(file.tellp()) + 3))
     258           0 :             n += 1024;
     259           7 :         n -= static_cast<index_t>(file.tellp()) + 3;
     260        6142 :         while (n > 0) {
     261        6135 :             file.put(' ');
     262        6135 :             n--;
     263        6135 :         }
     264             : 
     265             :         // close the header
     266           7 :         file << "\n}\n";
     267           7 :     }
     268             : 
     269             :     template <typename data_t>
     270             :     std::string EDF::getDataTypeName([[maybe_unused]] const DataContainer<data_t>& data)
     271           0 :     {
     272           0 :         throw InvalidArgumentError("EDF::getDataTypeName: invalid/unsupported data type");
     273           0 :     }
     274             : 
     275             :     template <>
     276             :     std::string EDF::getDataTypeName([[maybe_unused]] const DataContainer<int8_t>& data)
     277           0 :     {
     278           0 :         return "SignedByte";
     279           0 :     }
     280             :     template <>
     281             :     std::string EDF::getDataTypeName([[maybe_unused]] const DataContainer<uint8_t>& data)
     282           0 :     {
     283           0 :         return "UnsignedByte";
     284           0 :     }
     285             :     template <>
     286             :     std::string EDF::getDataTypeName([[maybe_unused]] const DataContainer<int16_t>& data)
     287           0 :     {
     288           0 :         return "SignedShort";
     289           0 :     }
     290             :     template <>
     291             :     std::string EDF::getDataTypeName([[maybe_unused]] const DataContainer<uint16_t>& data)
     292           0 :     {
     293           0 :         return "UnsignedShort";
     294           0 :     }
     295             :     template <>
     296             :     std::string EDF::getDataTypeName([[maybe_unused]] const DataContainer<float>& data)
     297           6 :     {
     298           6 :         return "FloatValue";
     299           6 :     }
     300             :     template <>
     301             :     std::string EDF::getDataTypeName([[maybe_unused]] const DataContainer<double>& data)
     302           1 :     {
     303           1 :         return "DoubleValue";
     304           1 :     }
     305             : 
     306             :     // ------------------------------------------
     307             :     // explicit template instantiation
     308             :     template DataContainer<float> EDF::read(std::string filename);
     309             :     template DataContainer<double> EDF::read(std::string filename);
     310             :     template DataContainer<index_t> EDF::read(std::string filename);
     311             : 
     312             :     template DataContainer<float> EDF::read(std::istream& input);
     313             :     template DataContainer<double> EDF::read(std::istream& input);
     314             :     template DataContainer<index_t> EDF::read(std::istream& input);
     315             : 
     316             :     template void EDF::write(const DataContainer<float>& data, std::string filename);
     317             :     template void EDF::write(const DataContainer<double>& data, std::string filename);
     318             :     template void EDF::write(const DataContainer<index_t>& data, std::string filename);
     319             : 
     320             :     template void EDF::write(const DataContainer<float>& data, std::ostream& output);
     321             :     template void EDF::write(const DataContainer<double>& data, std::ostream& output);
     322             :     template void EDF::write(const DataContainer<index_t>& data, std::ostream& output);
     323             : 
     324             :     template void EDF::writeHeader(std::ostream& file, const DataContainer<float>& data);
     325             :     template void EDF::writeHeader(std::ostream& file, const DataContainer<double>& data);
     326             :     template void EDF::writeHeader(std::ostream& file, const DataContainer<index_t>& data);
     327             : 
     328             : } // namespace elsa