allow user to set geometry size in z-direction

- add depth to geometry image dialog
- resolves #582

use/display compartment/membrane sizes in correct units

- use volume units for Compartment sizes
- use area units for Membrane sizes
- resolves #584

src/core/model/inc/model_compartments.hpp

@@ -26,6 +26,7 @@ class ModelGeometry;
 class ModelMembranes;
 class ModelSpecies;
 class ModelReactions;
+class ModelUnits;
 
 class ModelCompartments {
 private:
@@ -38,14 +39,15 @@ class ModelCompartments {
   ModelMembranes *modelMembranes = nullptr;
   ModelSpecies *modelSpecies = nullptr;
   ModelReactions *modelReactions = nullptr;
+  const ModelUnits *modelUnits{nullptr};
   simulate::SimulationData *simulationData = nullptr;
   bool hasUnsavedChanges{false};
 
 public:
   ModelCompartments();
   ModelCompartments(libsbml::Model *model, ModelGeometry *geometry,
                     ModelMembranes *membranes, ModelSpecies *species,
-                    ModelReactions *reactions, simulate::SimulationData *data);
+                    ModelReactions *reactions, const ModelUnits *units, simulate::SimulationData *data);
   const QStringList &getIds() const;
   const QStringList &getNames() const;
   const QVector<QRgb> &getColours() const;

src/core/model/inc/model_geometry.hpp

@@ -30,7 +30,9 @@ struct Settings;
 class ModelGeometry {
 private:
   double pixelWidth{1.0};
+  double pixelDepth{1.0};
   QPointF physicalOrigin{QPointF(0, 0)};
+  double zOrigin{0.0};
   QSizeF physicalSize{QSizeF(0, 0)};
   int numDimensions{3};
   QImage image;
@@ -45,6 +47,7 @@ class ModelGeometry {
   int importDimensions(const libsbml::Model *model);
   void convertSBMLGeometryTo3d();
   void writeDefaultGeometryToSBML();
+  void updateCompartmentAndMembraneSizes();
 
 public:
   ModelGeometry();
@@ -59,7 +62,10 @@ class ModelGeometry {
   void clear();
   int getNumDimensions() const;
   double getPixelWidth() const;
-  void setPixelWidth(double width);
+  void setPixelWidth(double width, bool updateSBML = true);
+  double getPixelDepth() const;
+  void setPixelDepth(double depth);
+  double getZOrigin() const;
   const QPointF &getPhysicalOrigin() const;
   const QSizeF &getPhysicalSize() const;
   const QImage &getImage() const;

src/core/model/inc/model_membranes.hpp

@@ -47,7 +47,7 @@ class ModelMembranes {
       const std::vector<std::unique_ptr<geometry::Compartment>> &compartments);
   void updateCompartmentImage(const QImage &img);
   void importMembraneIdsAndNames();
-  void exportToSBML(double pixelWidth);
+  void exportToSBML(double pixelArea);
   explicit ModelMembranes(libsbml::Model *model = nullptr);
   ModelMembranes(ModelMembranes &&) noexcept;
   ModelMembranes &operator=(ModelMembranes &&) noexcept;

src/core/model/inc/model_parameters.hpp

@@ -25,6 +25,7 @@ struct IdName {
 struct SpatialCoordinates {
   IdName x;
   IdName y;
+  IdName z;
 };
 
 struct IdNameValue {

src/core/model/src/model.cpp

@@ -6,6 +6,7 @@
 #include "utils.hpp"
 #include "validation.hpp"
 #include "xml_annotation.hpp"
+#include <QFileInfo>
 #include <algorithm>
 #include <sbml/SBMLTransforms.h>
 #include <sbml/SBMLTypes.h>
@@ -14,7 +15,6 @@
 #include <sbml/packages/spatial/extension/SpatialExtension.h>
 #include <stdexcept>
 #include <utility>
-#include <QFileInfo>
 
 namespace sme::model {
 
@@ -65,8 +65,9 @@ void Model::initModelData() {
   // matters, should be possible to reduce coupling here
   modelCompartments =
       ModelCompartments(model, &modelGeometry, &modelMembranes, &modelSpecies,
-                        &modelReactions, &getSimulationData());
-  modelGeometry = ModelGeometry(model, &modelCompartments, &modelMembranes, &settings);
+                        &modelReactions, &modelUnits, &getSimulationData());
+  modelGeometry =
+      ModelGeometry(model, &modelCompartments, &modelMembranes, &settings);
   modelGeometry.importSampledFieldGeometry(model);
   modelGeometry.importParametricGeometry(model, &settings);
   modelParameters = ModelParameters(model, &modelEvents);
@@ -151,7 +152,8 @@ void Model::exportSMEFile(const std::string &filename) {
 void Model::updateSBMLDoc() {
   modelGeometry.writeGeometryToSBML();
   setSbmlAnnotation(doc->getModel(), settings);
-  modelMembranes.exportToSBML(modelGeometry.getPixelWidth());
+  modelMembranes.exportToSBML(modelGeometry.getPixelWidth() *
+                              modelGeometry.getPixelDepth());
 }
 
 QString Model::getXml() {

src/core/model/src/model_compartments.cpp

@@ -6,6 +6,7 @@
 #include "model_membranes.hpp"
 #include "model_reactions.hpp"
 #include "model_species.hpp"
+#include "model_units.hpp"
 #include "sbml_utils.hpp"
 #include "simulate_data.hpp"
 #include <optional>
@@ -69,15 +70,14 @@ static void makeSizesValid(libsbml::Model *model) {
 
 ModelCompartments::ModelCompartments() = default;
 
-ModelCompartments::ModelCompartments(libsbml::Model *model,
-                                     ModelGeometry *geometry,
-                                     ModelMembranes *membranes,
-                                     ModelSpecies *species,
-                                     ModelReactions *reactions,
-                                     simulate::SimulationData *data)
+ModelCompartments::ModelCompartments(
+    libsbml::Model *model, ModelGeometry *geometry, ModelMembranes *membranes,
+    ModelSpecies *species, ModelReactions *reactions, const ModelUnits *units,
+    simulate::SimulationData *data)
     : ids{importIds(model)}, names{importNamesAndMakeUnique(model, ids)},
       sbmlModel{model}, modelGeometry{geometry}, modelMembranes{membranes},
-      modelSpecies{species}, modelReactions{reactions}, simulationData{data} {
+      modelSpecies{species}, modelReactions{reactions}, modelUnits{units},
+      simulationData{data} {
   makeSizesValid(model);
   colours = QVector<QRgb>(ids.size(), 0);
   compartments.reserve(static_cast<std::size_t>(ids.size()));
@@ -272,14 +272,17 @@ void ModelCompartments::setInteriorPoints(const QString &id,
   }
   const auto &origin{modelGeometry->getPhysicalOrigin()};
   auto pixelWidth{modelGeometry->getPixelWidth()};
+  auto centralZPoint{modelGeometry->getZOrigin() +
+                     modelGeometry->getPixelDepth() * 0.5};
   for (const auto &point : points) {
     SPDLOG_INFO("  - creating new interior point");
-    SPDLOG_INFO("    - pixel point: ({},{},{})", point.x(), point.y(), 0);
-    auto *ip = domain->createInteriorPoint();
+    SPDLOG_INFO("    - pixel point: ({},{},{})", point.x(), point.y(),
+                centralZPoint);
+    auto *ip{domain->createInteriorPoint()};
     // convert to physical units with pixelWidth and origin
     ip->setCoord1(origin.x() + pixelWidth * point.x());
     ip->setCoord2(origin.y() + pixelWidth * point.y());
-    ip->setCoord3(0.0);
+    ip->setCoord3(centralZPoint);
     SPDLOG_INFO("    - physical point: ({},{},0)", ip->getCoord1(),
                 ip->getCoord2(), ip->getCoord3());
   }
@@ -321,6 +324,10 @@ void ModelCompartments::setColour(const QString &id, QRgb colour) {
   geom->getDomainType(domainType)
       ->setSpatialDimensions(
           static_cast<int>(geom->getNumCoordinateComponents()));
+  if(compartment->isSetUnits()){
+    // we set the model units, compartment units are then inferred from that
+    compartment->unsetUnits();
+  }
   SPDLOG_INFO("  - sampledVolume '{}'", sfvol->getId());
   if (colour == 0 && sfvol->isSetSampledValue()) {
     sfvol->unsetSampledValue();
@@ -329,9 +336,16 @@ void ModelCompartments::setColour(const QString &id, QRgb colour) {
   }
   auto nPixels{compartments[static_cast<std::size_t>(i)]->nPixels()};
   double pixelWidth{modelGeometry->getPixelWidth()};
-  double volume{static_cast<double>(nPixels) * pixelWidth * pixelWidth};
-  compartment->setSize(volume);
-  SPDLOG_INFO("  - size {}", volume);
+  double pixelDepth{modelGeometry->getPixelDepth()};
+  double l3{static_cast<double>(nPixels) * pixelWidth * pixelWidth *
+            pixelDepth};
+  const auto &lengthUnit{modelUnits->getLength()};
+  const auto &volumeUnit{modelUnits->getVolume()};
+  double volOverL3 = getVolOverL3(lengthUnit, volumeUnit);
+  compartment->setSize(l3 / volOverL3);
+  SPDLOG_INFO("  - volume {} {}^3", l3, lengthUnit.name.toStdString());
+  SPDLOG_INFO("  - size {} {}", compartment->getSize(),
+              volumeUnit.name.toStdString());
   modelSpecies->updateCompartmentGeometry(id);
   modelMembranes->updateCompartments(compartments);
   modelMembranes->updateCompartmentNames(names);

src/core/model/src/model_compartments_t.cpp

@@ -40,13 +40,13 @@ SCENARIO("SBML compartments",
     REQUIRE(c.getIds().size() == 3);
     REQUIRE(c.getIds()[0] == "c1");
     REQUIRE(c.getName("c1") == "Outside");
-    REQUIRE(c.getSize("c1") == dbl_approx(5441));
+    REQUIRE(c.getSize("c1") == dbl_approx(5441000));
     REQUIRE(c.getIds()[1] == "c2");
     REQUIRE(c.getName("c2") == "Cell");
-    REQUIRE(c.getSize("c2") == dbl_approx(4034));
+    REQUIRE(c.getSize("c2") == dbl_approx(4034000));
     REQUIRE(c.getIds()[2] == "c3");
     REQUIRE(c.getName("c3") == "Nucleus");
-    REQUIRE(c.getSize("c3") == dbl_approx(525));
+    REQUIRE(c.getSize("c3") == dbl_approx(525000));
 
     // initial membranes
     REQUIRE(m.getIds().size() == 2);

src/core/model/src/model_geometry.cpp

@@ -28,13 +28,13 @@ int ModelGeometry::importDimensions(const libsbml::Model *model) {
       libsbml::CoordinateKind_t::SPATIAL_COORDINATEKIND_CARTESIAN_X)};
   if (xcoord == nullptr) {
     SPDLOG_ERROR("No x-coordinate found in SBML model");
-    return false;
+    return 0;
   }
   const auto *ycoord{geom->getCoordinateComponentByKind(
       libsbml::CoordinateKind_t::SPATIAL_COORDINATEKIND_CARTESIAN_Y)};
   if (ycoord == nullptr) {
     SPDLOG_ERROR("No y-coordinate found in SBML model");
-    return false;
+    return 0;
   }
 
   // import xy coordinates
@@ -48,19 +48,25 @@ int ModelGeometry::importDimensions(const libsbml::Model *model) {
   SPDLOG_INFO("  -> origin [{},{}]", physicalOrigin.x(), physicalOrigin.y());
   physicalSize = QSizeF(xmax - xmin, ymax - ymin);
   SPDLOG_INFO("  -> size [{},{}]", physicalSize.width(), physicalSize.height());
-
   if (nDim == 3) {
-    // todo: import z-direction min/max values
+    const auto *zcoord{geom->getCoordinateComponentByKind(
+        libsbml::CoordinateKind_t::SPATIAL_COORDINATEKIND_CARTESIAN_Z)};
+    if (zcoord == nullptr) {
+      SPDLOG_ERROR("No z-coordinate found in 3d SBML model");
+      return 0;
+    }
+    zOrigin = zcoord->getBoundaryMin()->getValue();
+    SPDLOG_INFO("  - found z origin {}", zOrigin);
+    pixelDepth = zcoord->getBoundaryMax()->getValue() - zOrigin;
+    SPDLOG_INFO("  - found z depth / pixel depth {}", pixelDepth);
   }
-
   return static_cast<int>(nDim);
 }
 
-static void createZCoordinateComponent(libsbml::Model *model){
-  auto *geom{
-      getOrCreateGeometry(model)};
-  auto *coordZ {geom->getCoordinateComponent(2)};
-  if(coordZ == nullptr){
+static void createZCoordinateComponent(libsbml::Model *model) {
+  auto *geom{getOrCreateGeometry(model)};
+  auto *coordZ{geom->getCoordinateComponent(2)};
+  if (coordZ == nullptr) {
     coordZ = geom->createCoordinateComponent();
   }
   coordZ->setType(
@@ -71,15 +77,15 @@ static void createZCoordinateComponent(libsbml::Model *model){
   paramZ->setUnits(model->getLengthUnits());
   paramZ->setConstant(true);
   paramZ->setValue(0);
-  auto *ssr {static_cast<libsbml::SpatialParameterPlugin *>(
-      paramZ->getPlugin("spatial"))
-      ->createSpatialSymbolReference()};
+  auto *ssr{static_cast<libsbml::SpatialParameterPlugin *>(
+                paramZ->getPlugin("spatial"))
+                ->createSpatialSymbolReference()};
   ssr->setSpatialRef(coordZ->getId());
   SPDLOG_INFO("  - creating Parameter: {}", paramZ->getId());
   SPDLOG_INFO("  - with spatialSymbolReference: {}", ssr->getSpatialRef());
   auto *minZ = coordZ->createBoundaryMin();
   minZ->setId("zBoundaryMin");
-  minZ->setValue(0);
+  minZ->setValue(0.0);
   auto *maxZ = coordZ->createBoundaryMax();
   maxZ->setId("zBoundaryMax");
   maxZ->setValue(1.0);
@@ -95,7 +101,8 @@ void ModelGeometry::convertSBMLGeometryTo3d() {
       return;
     }
     auto nDim{comp->getSpatialDimensions() + 1};
-    SPDLOG_INFO("Setting compartment '{}' dimensions to {}", comp->getId(), nDim);
+    SPDLOG_INFO("Setting compartment '{}' dimensions to {}", comp->getId(),
+                nDim);
     comp->setSpatialDimensions(nDim);
   }
 }
@@ -168,6 +175,15 @@ void ModelGeometry::writeDefaultGeometryToSBML() {
   createDefaultCompartmentGeometryIfMissing(sbmlModel);
 }
 
+void ModelGeometry::updateCompartmentAndMembraneSizes(){
+  // reassign all compartment colours to update sizes, interior points, etc
+  for (const auto &compartmentId : modelCompartments->getIds()) {
+    modelCompartments->setColour(compartmentId,
+                                 modelCompartments->getColour(compartmentId));
+  }
+  modelMembranes->exportToSBML(pixelWidth * pixelDepth);
+}
+
 ModelGeometry::ModelGeometry() = default;
 
 ModelGeometry::ModelGeometry(libsbml::Model *model,
@@ -179,7 +195,7 @@ ModelGeometry::ModelGeometry(libsbml::Model *model,
     SPDLOG_WARN("Failed to import geometry");
     writeDefaultGeometryToSBML();
     return;
-  } else if(nDim == 2){
+  } else if (nDim == 2) {
     convertSBMLGeometryTo3d();
     return;
   }
@@ -219,7 +235,7 @@ void ModelGeometry::importSampledFieldGeometry(const libsbml::Model *model) {
   image = gsf.image.convertToFormat(QImage::Format_Indexed8);
   hasImage = true;
   pixelWidth = calculatePixelWidth(image.size(), physicalSize);
-  setPixelWidth(pixelWidth);
+  setPixelWidth(pixelWidth, false);
   modelMembranes->updateCompartmentImage(image);
   for (const auto &[id, colour] : gsf.compartmentIdColourPairs) {
     SPDLOG_INFO("setting compartment {} colour to {:x}", id, colour);
@@ -334,7 +350,7 @@ int ModelGeometry::getNumDimensions() const { return numDimensions; }
 
 double ModelGeometry::getPixelWidth() const { return pixelWidth; }
 
-void ModelGeometry::setPixelWidth(double width) {
+void ModelGeometry::setPixelWidth(double width, bool updateSBML) {
   SPDLOG_INFO("Setting pixel width to {}", width);
   hasUnsavedChanges = true;
   double oldWidth = pixelWidth;
@@ -367,14 +383,31 @@ void ModelGeometry::setPixelWidth(double width) {
   min->setValue(physicalOrigin.y());
   max->setValue(physicalOrigin.y() + physicalSize.height());
   SPDLOG_INFO("  - y now in range [{},{}]", min->getValue(), max->getValue());
-
-  // reassign all compartment colours to update sizes, interior points, etc
-  for (const auto &compartmentId : modelCompartments->getIds()) {
-    modelCompartments->setColour(compartmentId,
-                                 modelCompartments->getColour(compartmentId));
+  if(updateSBML){
+    updateCompartmentAndMembraneSizes();
   }
 }
 
+double ModelGeometry::getPixelDepth() const { return pixelDepth; }
+
+void ModelGeometry::setPixelDepth(double depth) {
+  SPDLOG_INFO("Setting pixel depth to {}", depth);
+  pixelDepth = depth;
+  auto *geom{getOrCreateGeometry(sbmlModel)};
+  auto *coord{geom->getCoordinateComponentByKind(
+      libsbml::CoordinateKind_t::SPATIAL_COORDINATEKIND_CARTESIAN_Z)};
+  auto *min{coord->getBoundaryMin()};
+  auto *max{coord->getBoundaryMax()};
+  min->setValue(zOrigin);
+  max->setValue(zOrigin + pixelDepth);
+  SPDLOG_INFO("  - z now in range [{},{}]", min->getValue(), max->getValue());
+  updateCompartmentAndMembraneSizes();
+}
+
+double ModelGeometry::getZOrigin() const{
+  return zOrigin;
+}
+
 const QPointF &ModelGeometry::getPhysicalOrigin() const {
   return physicalOrigin;
 }

src/core/model/src/model_geometry_t.cpp

@@ -1,12 +1,7 @@
 #include "catch_wrapper.hpp"
 #include "mesh.hpp"
-#include "model_compartments.hpp"
+#include "model.hpp"
 #include "model_geometry.hpp"
-#include "model_membranes.hpp"
-#include "model_parameters.hpp"
-#include "model_reactions.hpp"
-#include "model_settings.hpp"
-#include "model_species.hpp"
 #include "serialization.hpp"
 #include <QFile>
 #include <sbml/SBMLTypes.h>
@@ -37,13 +32,14 @@ SCENARIO("Model geometry",
         libsbml::readSBMLFromString(f.readAll().toStdString().c_str()));
     model::Settings sbmlAnnotation;
     model::ModelCompartments mCompartments;
+    model::ModelUnits mUnits(doc->getModel());
     model::ModelMembranes mMembranes(doc->getModel());
     model::ModelGeometry mGeometry;
     model::ModelSpecies mSpecies;
     model::ModelReactions mReactions(doc->getModel(), &mMembranes);
     utils::SmeFileContents smeFileContents;
     mCompartments = model::ModelCompartments(
-        doc->getModel(), &mGeometry, &mMembranes, &mSpecies, &mReactions,
+        doc->getModel(), &mGeometry, &mMembranes, &mSpecies, &mReactions, &mUnits,
         &smeFileContents.simulationData);
     mGeometry = model::ModelGeometry(doc->getModel(), &mCompartments,
                                      &mMembranes, &sbmlAnnotation);