Merge pull request #92 from NeurodataWithoutBorders/add-spike-event-s…

…eries

Add SpikeEventSeries data type

CMakeLists.txt

@@ -27,6 +27,7 @@ add_library(
     src/nwb/base/TimeSeries.cpp
     src/nwb/device/Device.cpp
     src/nwb/ecephys/ElectricalSeries.cpp
+    src/nwb/ecephys/SpikeEventSeries.cpp
     src/nwb/file/ElectrodeGroup.cpp
     src/nwb/file/ElectrodeTable.cpp
     src/nwb/hdmf/base/Container.cpp

src/BaseIO.hpp

@@ -292,6 +292,14 @@ class BaseIO
   virtual std::unique_ptr<BaseRecordingData> getDataSet(
       const std::string& path) = 0;
 
+  /**
+   * @brief Checks whether a Dataset, Group, or Link already exists at the
+   * location in the file.
+   * @param path The location of the object in the file.
+   * @return Whether the object exists.
+   */
+  virtual bool objectExists(const std::string& path) = 0;
+
   /**
    * @brief Convenience function for creating NWB related attributes.
    * @param path The location of the object in the file.
@@ -325,6 +333,7 @@ class BaseIO
    * @return The status of the operation.
    */
   Status createTimestampsAttributes(const std::string& path);
+
   /**
    * @brief Returns true if the file is open.
    * @return True if the file is open, false otherwise.

src/hdf5/HDF5IO.cpp

@@ -437,6 +437,16 @@ bool HDF5IO::canModifyObjects()
   return statusOK && !inSWMRMode;
 }
 
+bool HDF5IO::objectExists(const std::string& path)
+{
+  htri_t exists = H5Lexists(file->getId(), path.c_str(), H5P_DEFAULT);
+  if (exists > 0) {
+    return true;
+  } else {
+    return false;
+  }
+}
+
 std::unique_ptr<AQNWB::BaseRecordingData> HDF5IO::getDataSet(
     const std::string& path)
 {

src/hdf5/HDF5IO.hpp

@@ -234,6 +234,14 @@ class HDF5IO : public BaseIO
   std::unique_ptr<BaseRecordingData> getDataSet(
       const std::string& path) override;
 
+  /**
+   * @brief Checks whether a Dataset, Group, or Link already exists at the
+   * location in the file.
+   * @param path The location of the object in the file.
+   * @return Whether the object exists.
+   */
+  bool objectExists(const std::string& path) override;
+
   /**
    * @brief Returns the HDF5 type of object at a given path.
    * @param path The location in the file of the object.

src/nwb/NWBFile.cpp

@@ -13,15 +13,18 @@
 #include "Utils.hpp"
 #include "nwb/device/Device.hpp"
 #include "nwb/ecephys/ElectricalSeries.hpp"
+#include "nwb/ecephys/SpikeEventSeries.hpp"
 #include "nwb/file/ElectrodeGroup.hpp"
-#include "nwb/file/ElectrodeTable.hpp"
 #include "spec/core.hpp"
 #include "spec/hdmf_common.hpp"
 #include "spec/hdmf_experimental.hpp"
 
 using namespace AQNWB::NWB;
 
-constexpr SizeType CHUNK_XSIZE = 2048;
+constexpr SizeType CHUNK_XSIZE =
+    2048;  // TODO - replace these with io settings input
+constexpr SizeType SPIKE_CHUNK_XSIZE =
+    8;  // TODO - replace with io settings input
 
 std::vector<SizeType> NWBFile::emptyContainerIndexes = {};
 
@@ -99,6 +102,7 @@ Status NWBFile::createFileStructure(std::string description,
 
 Status NWBFile::createElectricalSeries(
     std::vector<Types::ChannelVector> recordingArrays,
+    std::vector<std::string> recordingNames,
     const BaseDataType& dataType,
     RecordingContainers* recordingContainers,
     std::vector<SizeType>& containerIndexes)
@@ -107,27 +111,44 @@ Status NWBFile::createElectricalSeries(
     return Status::Failure;
   }
 
-  // store all recorded data in the acquisition group
-  std::string rootPath = "/acquisition/";
+  if (recordingNames.size() != recordingArrays.size()) {
+    return Status::Failure;
+  }
+
+  // Setup electrode table if it was not yet created
+  bool electrodeTableCreated =
+      io->objectExists(ElectrodeTable::electrodeTablePath);
+  if (!electrodeTableCreated) {
+    elecTable = std::make_unique<ElectrodeTable>(io);
+    elecTable->initialize();
 
-  // Setup electrode table
-  ElectrodeTable elecTable = ElectrodeTable(io);
-  elecTable.initialize();
+    // Add electrode information to table (does not write to datasets yet)
+    for (const auto& channelVector : recordingArrays) {
+      elecTable->addElectrodes(channelVector);
+    }
+  }
+
+  // Create datasets
+  for (size_t i = 0; i < recordingArrays.size(); ++i) {
+    const auto& channelVector = recordingArrays[i];
+    const std::string& recordingName = recordingNames[i];
 
-  // Create continuous datasets
-  for (const auto& channelVector : recordingArrays) {
     // Setup electrodes and devices
     std::string groupName = channelVector[0].groupName;
     std::string devicePath = "/general/devices/" + groupName;
     std::string electrodePath = "/general/extracellular_ephys/" + groupName;
-    std::string electricalSeriesPath = rootPath + groupName;
+    std::string electricalSeriesPath = acquisitionPath + "/" + recordingName;
 
-    Device device = Device(devicePath, io, "description", "unknown");
-    device.initialize();
+    // Check if device exists for groupName, create device and electrode group
+    // if not
+    if (!io->objectExists(devicePath)) {
+      Device device = Device(devicePath, io, "description", "unknown");
+      device.initialize();
 
-    ElectrodeGroup elecGroup =
-        ElectrodeGroup(electrodePath, io, "description", "unknown", device);
-    elecGroup.initialize();
+      ElectrodeGroup elecGroup =
+          ElectrodeGroup(electrodePath, io, "description", "unknown", device);
+      elecGroup.initialize();
+    }
 
     // Setup electrical series datasets
     auto electricalSeries = std::make_unique<ElectricalSeries>(
@@ -142,14 +163,96 @@ Status NWBFile::createElectricalSeries(
     electricalSeries->initialize();
     recordingContainers->addContainer(std::move(electricalSeries));
     containerIndexes.push_back(recordingContainers->containers.size() - 1);
+  }
 
-    // Add electrode information to electrode table (does not write to datasets
-    // yet)
-    elecTable.addElectrodes(channelVector);
+  // write electrode information to datasets
+  // (requires that the ElectrodeGroup has been written)
+  if (!electrodeTableCreated) {
+    elecTable->finalize();
+  }
+
+  return Status::Success;
+}
+
+Status NWBFile::createSpikeEventSeries(
+    std::vector<Types::ChannelVector> recordingArrays,
+    std::vector<std::string> recordingNames,
+    const BaseDataType& dataType,
+    RecordingContainers* recordingContainers,
+    std::vector<SizeType>& containerIndexes)
+{
+  if (!io->canModifyObjects()) {
+    return Status::Failure;
+  }
+
+  if (recordingNames.size() != recordingArrays.size()) {
+    return Status::Failure;
+  }
+
+  // Setup electrode table if it was not yet created
+  bool electrodeTableCreated =
+      io->objectExists(ElectrodeTable::electrodeTablePath);
+  if (!electrodeTableCreated) {
+    elecTable = std::make_unique<ElectrodeTable>(io);
+    elecTable->initialize();
+
+    // Add electrode information to table (does not write to datasets yet)
+    for (const auto& channelVector : recordingArrays) {
+      elecTable->addElectrodes(channelVector);
+    }
+  }
+
+  // Create datasets
+  for (size_t i = 0; i < recordingArrays.size(); ++i) {
+    const auto& channelVector = recordingArrays[i];
+    const std::string& recordingName = recordingNames[i];
+
+    // Setup electrodes and devices
+    std::string groupName = channelVector[0].groupName;
+    std::string devicePath = "/general/devices/" + groupName;
+    std::string electrodePath = "/general/extracellular_ephys/" + groupName;
+    std::string spikeEventSeriesPath = acquisitionPath + "/" + recordingName;
+
+    // Check if device exists for groupName, create device and electrode group
+    // if not
+    if (!io->objectExists(devicePath)) {
+      Device device = Device(devicePath, io, "description", "unknown");
+      device.initialize();
+
+      ElectrodeGroup elecGroup =
+          ElectrodeGroup(electrodePath, io, "description", "unknown", device);
+      elecGroup.initialize();
+    }
+
+    // Setup Spike Event Series datasets
+    SizeArray dsetSize;
+    SizeArray chunkSize;
+    if (channelVector.size() == 1) {
+      dsetSize = SizeArray {0, 0};
+      chunkSize = SizeArray {SPIKE_CHUNK_XSIZE, 1};
+    } else {
+      dsetSize = SizeArray {0, channelVector.size(), 0};
+      chunkSize = SizeArray {SPIKE_CHUNK_XSIZE, 1, 1};
+    }
+
+    auto spikeEventSeries = std::make_unique<SpikeEventSeries>(
+        spikeEventSeriesPath,
+        io,
+        dataType,
+        channelVector,
+        "Stores spike waveforms from an extracellular ephys recording",
+        dsetSize,
+        chunkSize);
+    spikeEventSeries->initialize();
+    recordingContainers->addContainer(std::move(spikeEventSeries));
+    containerIndexes.push_back(recordingContainers->containers.size() - 1);
   }
 
   // write electrode information to datasets
-  elecTable.finalize();
+  // (requires that the ElectrodeGroup has been written)
+  if (!electrodeTableCreated) {
+    elecTable->finalize();
+  }
 
   return Status::Success;
 }
@@ -161,7 +264,7 @@ void NWBFile::cacheSpecifications(
     const std::array<std::pair<std::string_view, std::string_view>, N>&
         specVariables)
 {
-  io->createGroup("/specifications/" + specPath + "/");
+  io->createGroup("/specifications/" + specPath);
   io->createGroup("/specifications/" + specPath + "/" + versionNumber);
 
   for (const auto& [name, content] : specVariables) {

src/nwb/NWBFile.hpp

@@ -11,6 +11,7 @@
 #include "Types.hpp"
 #include "nwb/RecordingContainers.hpp"
 #include "nwb/base/TimeSeries.hpp"
+#include "nwb/file/ElectrodeTable.hpp"
 
 /*!
  * \namespace AQNWB::NWB
@@ -71,14 +72,38 @@ class NWBFile
    * @param recordingArrays vector of ChannelVector indicating the electrodes to
    *                        record from. A separate ElectricalSeries will be
    *                        created for each ChannelVector.
+   * @param recordingNames vector indicating the names of the ElectricalSeries
+   * within the acquisition group
+   * @param dataType The data type of the elements in the data block.
    * @param recordingContainers The container to store the created TimeSeries.
    * @param containerIndexes The indexes of the containers added to
    * recordingContainers
-   * @param dataType The data type of the elements in the data block.
    * @return Status The status of the object creation operation.
    */
   Status createElectricalSeries(
       std::vector<Types::ChannelVector> recordingArrays,
+      std::vector<std::string> recordingNames,
+      const BaseDataType& dataType = BaseDataType::I16,
+      RecordingContainers* recordingContainers = nullptr,
+      std::vector<SizeType>& containerIndexes = emptyContainerIndexes);
+
+  /**
+   * @brief Create SpikeEventSeries objects to record data into.
+   * Created objects are stored in recordingContainers.
+   * @param recordingArrays vector of ChannelVector indicating the electrodes to
+   *                        record from. A separate ElectricalSeries will be
+   *                        created for each ChannelVector.
+   * @param recordingNames vector indicating the names of the SpikeEventSeries
+   * within the acquisition group
+   * @param dataType The data type of the elements in the data block.
+   * @param recordingContainers The container to store the created TimeSeries.
+   * @param containerIndexes The indexes of the containers added to
+   * recordingContainers
+   * @return Status The status of the object creation operation.
+   */
+  Status createSpikeEventSeries(
+      std::vector<Types::ChannelVector> recordingArrays,
+      std::vector<std::string> recordingNames,
       const BaseDataType& dataType = BaseDataType::I16,
       RecordingContainers* recordingContainers = nullptr,
       std::vector<SizeType>& containerIndexes = emptyContainerIndexes);
@@ -127,9 +152,11 @@ class NWBFile
       const std::array<std::pair<std::string_view, std::string_view>, N>&
           specVariables);
 
+  std::unique_ptr<ElectrodeTable> elecTable;
   const std::string identifierText;
   std::shared_ptr<BaseIO> io;
   static std::vector<SizeType> emptyContainerIndexes;
+  inline const static std::string acquisitionPath = "/acquisition";
 };
 
 }  // namespace AQNWB::NWB

src/nwb/RecordingContainers.cpp

@@ -2,6 +2,7 @@
 #include "nwb/RecordingContainers.hpp"
 
 #include "nwb/ecephys/ElectricalSeries.hpp"
+#include "nwb/ecephys/SpikeEventSeries.hpp"
 #include "nwb/hdmf/base/Container.hpp"
 
 using namespace AQNWB::NWB;
@@ -63,3 +64,18 @@ Status RecordingContainers::writeElectricalSeriesData(
 
   es->writeChannel(channel.localIndex, numSamples, data, timestamps);
 }
+
+Status RecordingContainers::writeSpikeEventData(const SizeType& containerInd,
+                                                const SizeType& numSamples,
+                                                const SizeType& numChannels,
+                                                const void* data,
+                                                const void* timestamps)
+{
+  SpikeEventSeries* ses =
+      dynamic_cast<SpikeEventSeries*>(getContainer(containerInd));
+
+  if (ses == nullptr)
+    return Status::Failure;
+
+  ses->writeSpike(numSamples, numChannels, data, timestamps);
+}

src/nwb/RecordingContainers.hpp

@@ -71,7 +71,7 @@ class RecordingContainers
                              const void* timestamps);
 
   /**
-   * @brief Write ElectricalSereis data to a recordingContainer dataset.
+   * @brief Write ElectricalSeries data to a recordingContainer dataset.
    * @param containerInd The index of the electrical series dataset within the
    * electrical series group.
    * @param channel The channel index to use for writing timestamps.
@@ -89,6 +89,22 @@ class RecordingContainers
                                    const void* data,
                                    const void* timestamps);
 
+  /**
+   * @brief Write SpikeEventSeries data to a recordingContainer dataset.
+   * @param containerInd The index of the SpikeEventSeries dataset within the
+   * SpikeEventSeries containers.
+   * @param numSamples Number of samples in the time for the single event.
+   * @param numChannels Number of channels in the time for the single event.
+   * @param data A pointer to the data block.
+   * @param timestamps A pointer to the timestamps block
+   * @return The status of the write operation.
+   */
+  Status writeSpikeEventData(const SizeType& containerInd,
+                             const SizeType& numSamples,
+                             const SizeType& numChannels,
+                             const void* data,
+                             const void* timestamps);
+
   std::vector<std::unique_ptr<Container>> containers;
   std::string name;
 };