InMemorySerializer with endianess swap.
authorChristian Blau <cblau@gwdg.de>
Tue, 29 Oct 2019 14:42:46 +0000 (15:42 +0100)
committerMark Abraham <mark.j.abraham@gmail.com>
Sun, 24 Nov 2019 22:00:31 +0000 (23:00 +0100)
Adds an InMemorySerializer that swaps endianess. Required for reading
data that is not handled with XDR serializers and was written with
different endianess, e.g., mrc files for density guided simulations.

refs #2282

Change-Id: Ice96fe8684bb0b920c22dc7d3c2b2087260eb9ad

src/gromacs/utility/inmemoryserializer.cpp
src/gromacs/utility/inmemoryserializer.h
src/gromacs/utility/tests/CMakeLists.txt
src/gromacs/utility/tests/inmemoryserializer.cpp [new file with mode: 0644]

index c1adcee5eefb37ad837ae572f4565fd3f2f2537d..eb012c36fccb4a9fe44d3fa47a848ed5510a72ab 100644 (file)
@@ -39,8 +39,6 @@
 #include <algorithm>
 #include <vector>
 
-#include "gromacs/utility/gmxassert.h"
-
 namespace gmx
 {
 
@@ -51,7 +49,7 @@ template<typename T>
 class CharBuffer
 {
 public:
-    static const size_t ValueSize = sizeof(T);
+    static constexpr size_t ValueSize = sizeof(T);
 
     explicit CharBuffer(T value) { u.v = value; }
     explicit CharBuffer(const char buffer[]) { std::copy(buffer, buffer + ValueSize, u.c); }
@@ -70,6 +68,26 @@ private:
     } u;
 };
 
+//! Return \c value with the byte order swapped.
+template<typename T>
+T swapEndian(const T& value)
+{
+    union {
+        T                           value_;
+        std::array<char, sizeof(T)> valueAsCharArray_;
+    } endianessSwappedValue;
+
+    endianessSwappedValue.value_ = value;
+    int hiByte                   = sizeof(T) - 1;
+    for (int loByte = 0; hiByte > loByte; loByte++, hiByte--)
+    {
+        std::swap(endianessSwappedValue.valueAsCharArray_[loByte],
+                  endianessSwappedValue.valueAsCharArray_[hiByte]);
+    }
+
+    return endianessSwappedValue.value_;
+}
+
 } // namespace
 
 /********************************************************************
@@ -79,10 +97,18 @@ private:
 class InMemorySerializer::Impl
 {
 public:
+    Impl(EndianSwapBehavior endianSwapBehavior) : endianSwapBehavior_(endianSwapBehavior) {}
     template<typename T>
     void doValue(T value)
     {
-        CharBuffer<T>(value).appendTo(&buffer_);
+        if (endianSwapBehavior_ == EndianSwapBehavior::DoSwap)
+        {
+            CharBuffer<T>(swapEndian(value)).appendTo(&buffer_);
+        }
+        else
+        {
+            CharBuffer<T>(value).appendTo(&buffer_);
+        }
     }
     void doString(const std::string& value)
     {
@@ -90,12 +116,16 @@ public:
         buffer_.insert(buffer_.end(), value.begin(), value.end());
     }
 
-    std::vector<char> buffer_;
+    std::vector<char>  buffer_;
+    EndianSwapBehavior endianSwapBehavior_;
 };
 
-InMemorySerializer::InMemorySerializer() : impl_(new Impl) {}
+InMemorySerializer::InMemorySerializer(EndianSwapBehavior endianSwapBehavior) :
+    impl_(new Impl(endianSwapBehavior))
+{
+}
 
-InMemorySerializer::~InMemorySerializer() {}
+InMemorySerializer::~InMemorySerializer() = default;
 
 std::vector<char> InMemorySerializer::finishAndGetBuffer()
 {
@@ -180,17 +210,25 @@ void InMemorySerializer::doString(std::string* value)
 class InMemoryDeserializer::Impl
 {
 public:
-    explicit Impl(ArrayRef<const char> buffer, bool sourceIsDouble) :
+    explicit Impl(ArrayRef<const char> buffer, bool sourceIsDouble, EndianSwapBehavior endianSwapBehavior) :
         buffer_(buffer),
         sourceIsDouble_(sourceIsDouble),
-        pos_(0)
+        pos_(0),
+        endianSwapBehavior_(endianSwapBehavior)
     {
     }
 
     template<typename T>
     void doValue(T* value)
     {
-        *value = CharBuffer<T>(&buffer_[pos_]).value();
+        if (endianSwapBehavior_ == EndianSwapBehavior::DoSwap)
+        {
+            *value = swapEndian(CharBuffer<T>(&buffer_[pos_]).value());
+        }
+        else
+        {
+            *value = CharBuffer<T>(&buffer_[pos_]).value();
+        }
         pos_ += CharBuffer<T>::ValueSize;
     }
     void doString(std::string* value)
@@ -204,14 +242,17 @@ public:
     ArrayRef<const char> buffer_;
     bool                 sourceIsDouble_;
     size_t               pos_;
+    EndianSwapBehavior   endianSwapBehavior_;
 };
 
-InMemoryDeserializer::InMemoryDeserializer(ArrayRef<const char> buffer, bool sourceIsDouble) :
-    impl_(new Impl(buffer, sourceIsDouble))
+InMemoryDeserializer::InMemoryDeserializer(ArrayRef<const char> buffer,
+                                           bool                 sourceIsDouble,
+                                           EndianSwapBehavior   endianSwapBehavior) :
+    impl_(new Impl(buffer, sourceIsDouble, endianSwapBehavior))
 {
 }
 
-InMemoryDeserializer::~InMemoryDeserializer() {}
+InMemoryDeserializer::~InMemoryDeserializer() = default;
 
 bool InMemoryDeserializer::sourceIsDouble() const
 {
index 60998c15aa75ac7b93e8dbbd8bb8eeb29e75a77d..8b6ac98a70f1ba35430753dc3f2fe3ada0cf4afa 100644 (file)
 namespace gmx
 {
 
+//! Specify endian swapping behvaior
+enum class EndianSwapBehavior : int
+{
+    DoSwap,    //!< Swap the bytes
+    DoNotSwap, //!< Do not swap the bytes
+    Count      //!< Number of possible behaviors
+};
+
 class InMemorySerializer : public ISerializer
 {
 public:
-    InMemorySerializer();
+    explicit InMemorySerializer(EndianSwapBehavior endianSwapBehavior = EndianSwapBehavior::DoNotSwap);
     ~InMemorySerializer() override;
 
     std::vector<char> finishAndGetBuffer();
@@ -85,7 +93,9 @@ private:
 class InMemoryDeserializer : public ISerializer
 {
 public:
-    explicit InMemoryDeserializer(ArrayRef<const char> buffer, bool sourceIsDouble);
+    InMemoryDeserializer(ArrayRef<const char> buffer,
+                         bool                 sourceIsDouble,
+                         EndianSwapBehavior   endianSwapBehavior = EndianSwapBehavior::DoNotSwap);
     ~InMemoryDeserializer() override;
 
     //! Get if the source data was written in double precsion
index ea9bfe398ba358dfaab1284f282c2bfc0739d28f..e8040f05f7a5909654efe9add4ab6402d3b91263 100644 (file)
@@ -40,6 +40,7 @@ gmx_add_unit_test(UtilityUnitTests utility-test
                   defaultinitializationallocator.cpp
                   enumerationhelpers.cpp
                   fixedcapacityvector.cpp
+                  inmemoryserializer.cpp
                   keyvaluetreeserializer.cpp
                   keyvaluetreetransform.cpp
                   logger.cpp
diff --git a/src/gromacs/utility/tests/inmemoryserializer.cpp b/src/gromacs/utility/tests/inmemoryserializer.cpp
new file mode 100644 (file)
index 0000000..a3612eb
--- /dev/null
@@ -0,0 +1,246 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2019, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
+ * and including many others, as listed in the AUTHORS file in the
+ * top-level source directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+/*! \internal \file
+ * \brief
+ * Tests for gmx::InMemorySerializer and InMemoryDeserializer.
+ *
+ * \author Christian Blau <blau@kth.se>
+ * \ingroup module_utility
+ */
+
+#include "gmxpre.h"
+
+#include "gromacs/utility/inmemoryserializer.h"
+
+#include <gmock/gmock.h>
+#include <gtest/gtest.h>
+
+namespace gmx
+{
+namespace test
+{
+namespace
+{
+union IntAndFloat32 {
+    std::int32_t int32Value_;
+    float        floatValue_;
+};
+
+union IntAndFloat64 {
+    std::int64_t int64Value_;
+    double       doubleValue_;
+};
+
+//! Constants used for testing endian swap operations
+//! \{
+constexpr std::int16_t c_int16Value        = static_cast<std::int16_t>(0x7A2B);
+constexpr std::int16_t c_int16ValueSwapped = static_cast<std::int16_t>(0x2B7A);
+constexpr std::int32_t c_int32Value        = static_cast<std::int32_t>(0x78ABCDEF);
+constexpr std::int32_t c_int32ValueSwapped = static_cast<std::int32_t>(0xEFCDAB78);
+constexpr std::int64_t c_int64Value        = static_cast<std::int64_t>(0x78ABCDEF12345678);
+constexpr std::int64_t c_int64ValueSwapped = static_cast<std::int64_t>(0x78563412EFCDAB78);
+
+constexpr const IntAndFloat32 c_intAndFloat32{ c_int32Value };
+constexpr const IntAndFloat64 c_intAndFloat64{ c_int64Value };
+
+constexpr const IntAndFloat32 c_intAndFloat32Swapped{ c_int32ValueSwapped };
+constexpr const IntAndFloat64 c_intAndFloat64Swapped{ c_int64ValueSwapped };
+//! \}
+
+//! Return the integer used for testing, depending on the size of int.
+constexpr int integerSizeDependentTestingValue()
+{
+    return sizeof(int) == 4 ? c_int32Value : sizeof(int) == 8 ? c_int64Value : c_int16Value;
+}
+
+//! Return the endianess-swapped integer used for testing, depending on the size of int.
+constexpr int integerSizeDependentTestingValueEndianessSwapped()
+{
+    return sizeof(int) == 4 ? c_int32ValueSwapped
+                            : sizeof(int) == 8 ? c_int64ValueSwapped : c_int16ValueSwapped;
+}
+
+class InMemorySerializerTest : public ::testing::Test
+{
+public:
+    struct SerializerValues
+    {
+        bool           boolValue_;
+        unsigned char  unsignedCharValue_;
+        char           charValue_;
+        unsigned short unsignedShortValue_;
+        std::int32_t   int32Value_;
+        float          floatValue_;
+        std::int64_t   int64Value_;
+        double         doubleValue_;
+        int            intValue_;
+        real           realValue_;
+    };
+
+    void serialize(ISerializer* serializer, SerializerValues* values)
+    {
+        EXPECT_FALSE(serializer->reading());
+        doValues(serializer, values);
+    }
+
+    SerializerValues deserialize(ISerializer* serializer)
+    {
+        EXPECT_TRUE(serializer->reading());
+        SerializerValues result;
+        doValues(serializer, &result);
+        return result;
+    }
+
+    void checkSerializerValuesforEquality(const SerializerValues& lhs, const SerializerValues& rhs)
+    {
+        EXPECT_EQ(lhs.boolValue_, rhs.boolValue_);
+        EXPECT_EQ(lhs.unsignedCharValue_, rhs.unsignedCharValue_);
+        EXPECT_EQ(lhs.charValue_, rhs.charValue_);
+        EXPECT_EQ(lhs.intValue_, rhs.intValue_);
+        EXPECT_EQ(lhs.int32Value_, rhs.int32Value_);
+        EXPECT_EQ(lhs.int64Value_, rhs.int64Value_);
+        EXPECT_EQ(lhs.unsignedShortValue_, rhs.unsignedShortValue_);
+        EXPECT_EQ(lhs.realValue_, rhs.realValue_);
+        EXPECT_EQ(lhs.floatValue_, rhs.floatValue_);
+        EXPECT_EQ(lhs.doubleValue_, rhs.doubleValue_);
+    }
+
+private:
+    void doValues(ISerializer* serializer, SerializerValues* values)
+    {
+        serializer->doBool(&values->boolValue_);
+        serializer->doUChar(&values->unsignedCharValue_);
+        serializer->doChar(&values->charValue_);
+        serializer->doInt(&values->intValue_);
+        serializer->doInt32(&values->int32Value_);
+        serializer->doInt64(&values->int64Value_);
+        serializer->doUShort(&values->unsignedShortValue_);
+        serializer->doReal(&values->realValue_);
+        serializer->doFloat(&values->floatValue_);
+        serializer->doDouble(&values->doubleValue_);
+    }
+
+protected:
+    SerializerValues defaultValues_ = { true,
+                                        0x78,
+                                        0x78,
+                                        static_cast<unsigned short>(c_int16Value),
+                                        c_int32Value,
+                                        c_intAndFloat32.floatValue_,
+                                        c_int64Value,
+                                        c_intAndFloat64.doubleValue_,
+                                        integerSizeDependentTestingValue(),
+                                        std::is_same<real, double>::value
+                                                ? static_cast<real>(c_intAndFloat64.doubleValue_)
+                                                : static_cast<real>(c_intAndFloat32.floatValue_) };
+
+    SerializerValues endianessSwappedValues_ = {
+        true,
+        0x78,
+        0x78,
+        static_cast<unsigned short>(c_int16ValueSwapped),
+        c_int32ValueSwapped,
+        c_intAndFloat32Swapped.floatValue_,
+        c_int64ValueSwapped,
+        c_intAndFloat64Swapped.doubleValue_,
+        integerSizeDependentTestingValueEndianessSwapped(),
+        std::is_same<real, float>::value ? static_cast<real>(c_intAndFloat32Swapped.floatValue_)
+                                         : static_cast<real>(c_intAndFloat64Swapped.doubleValue_)
+    };
+};
+
+TEST_F(InMemorySerializerTest, Roundtrip)
+{
+    InMemorySerializer serializer;
+    SerializerValues   values = defaultValues_;
+    serialize(&serializer, &values);
+
+    auto buffer = serializer.finishAndGetBuffer();
+
+    InMemoryDeserializer deserializer(buffer, std::is_same<real, double>::value);
+
+    SerializerValues deserialisedValues = deserialize(&deserializer);
+
+    checkSerializerValuesforEquality(values, deserialisedValues);
+}
+
+TEST_F(InMemorySerializerTest, RoundtripWithEndianessSwap)
+{
+    InMemorySerializer serializerWithSwap(EndianSwapBehavior::DoSwap);
+    SerializerValues   values = defaultValues_;
+    serialize(&serializerWithSwap, &values);
+
+    auto buffer = serializerWithSwap.finishAndGetBuffer();
+
+    InMemoryDeserializer deserializerWithSwap(buffer, std::is_same<real, double>::value,
+                                              EndianSwapBehavior::DoSwap);
+
+    SerializerValues deserialisedValues = deserialize(&deserializerWithSwap);
+
+    checkSerializerValuesforEquality(values, deserialisedValues);
+}
+
+TEST_F(InMemorySerializerTest, SerializerExplicitEndianessSwap)
+{
+    InMemorySerializer serializerWithSwap(EndianSwapBehavior::DoSwap);
+    SerializerValues   values = defaultValues_;
+    serialize(&serializerWithSwap, &values);
+
+    auto buffer = serializerWithSwap.finishAndGetBuffer();
+
+    InMemoryDeserializer deserializerWithOutSwap(buffer, std::is_same<real, double>::value);
+
+    SerializerValues deserialisedValues = deserialize(&deserializerWithOutSwap);
+    checkSerializerValuesforEquality(endianessSwappedValues_, deserialisedValues);
+}
+
+TEST_F(InMemorySerializerTest, DeserializerExplicitEndianessSwap)
+{
+    InMemorySerializer serializer;
+    SerializerValues   values = defaultValues_;
+    serialize(&serializer, &values);
+
+    auto buffer = serializer.finishAndGetBuffer();
+
+    InMemoryDeserializer deserializerWithSwap(buffer, std::is_same<real, double>::value,
+                                              EndianSwapBehavior::DoSwap);
+
+    SerializerValues deserialisedValues = deserialize(&deserializerWithSwap);
+    checkSerializerValuesforEquality(endianessSwappedValues_, deserialisedValues);
+}
+
+} // namespace
+} // namespace test
+} // namespace gmx