Skip to content

Commit

Permalink
Merge pull request #1202 from virtualcell/FixingUnits
Browse files Browse the repository at this point in the history
Fixed Issue where Geometry did not scale with Unit change to length
  • Loading branch information
CodeByDrescher committed Apr 1, 2024
2 parents 238a954 + 39284b0 commit 4c0c15b
Show file tree
Hide file tree
Showing 23 changed files with 282 additions and 247 deletions.
Expand Up @@ -16,6 +16,7 @@
import javax.swing.UIManager;
import javax.swing.UnsupportedLookAndFeelException;

import cbit.vcell.resource.ResourceUtil;
import org.apache.logging.log4j.LogManager;
import org.apache.logging.log4j.Logger;

Expand All @@ -32,7 +33,9 @@ public static void setVCellLookAndFeel() {
try {
lg.info("Operating system: " + osi.getOsType());
lg.info("About to set the look and feel. Before setting, we're using: " + UIManager.getLookAndFeel().getName());
UIManager.setLookAndFeel(UIManager.getSystemLookAndFeelClassName());
String lookAndFeelType = OperatingSystemInfo.getInstance().isLinux() ?
UIManager.getCrossPlatformLookAndFeelClassName() : UIManager.getSystemLookAndFeelClassName();
UIManager.setLookAndFeel(lookAndFeelType);
} catch (ClassNotFoundException | InstantiationException | IllegalAccessException | UnsupportedLookAndFeelException e) {
lg.warn("Error while setting look and feel:", e);
}
Expand Down
131 changes: 101 additions & 30 deletions vcell-core/src/main/java/cbit/vcell/biomodel/ModelUnitConverter.java
@@ -1,20 +1,26 @@
package cbit.vcell.biomodel;

import cbit.image.ImageException;
import cbit.vcell.geometry.*;
import cbit.vcell.geometry.surface.GeometrySurfaceDescription;
import cbit.vcell.mapping.*;
import cbit.vcell.mapping.BioEvent.EventAssignment;
import cbit.vcell.math.MathDescription;
import cbit.vcell.matrix.RationalExp;
import cbit.vcell.matrix.RationalNumber;
import cbit.vcell.model.*;
import cbit.vcell.parser.*;
import cbit.vcell.solver.Simulation;
import cbit.vcell.units.VCUnitDefinition;
import cbit.vcell.xml.XMLSource;
import cbit.vcell.xml.XmlHelper;
import cbit.vcell.xml.XmlParseException;
import org.apache.logging.log4j.LogManager;
import org.apache.logging.log4j.Logger;
import org.vcell.util.Extent;
import org.vcell.util.Origin;

import java.math.BigDecimal;
import java.math.RoundingMode;
import java.util.LinkedHashMap;
import java.util.Map;

Expand All @@ -24,25 +30,25 @@ public class ModelUnitConverter {

public static ModelUnitSystem createSbmlModelUnitSystem() {
final String substanceUnit = "umol";
String volumeSubstanceSymbol = substanceUnit;
String membraneSubstanceSymbol = substanceUnit;
String lumpedReactionSubstanceSymbol = substanceUnit;
String volumeSymbol = "l";
// All three of the following are the same in an SBML Unit system
String volumeSubstanceSymbol = substanceUnit; // Used to demonstrate difference between VCell and SBML
String membraneSubstanceSymbol = substanceUnit; // Used to demonstrate difference between VCell and SBML
String lumpedReactionSubstanceSymbol = substanceUnit; // Used to demonstrate difference between VCell and SBML
String volumeSymbol = "l";
String areaSymbol = "dm2";
String lengthSymbol = "dm";
String timeSymbol = "s";
ModelUnitSystem mus = ModelUnitSystem.createVCModelUnitSystem(volumeSubstanceSymbol, membraneSubstanceSymbol, lumpedReactionSubstanceSymbol, volumeSymbol, areaSymbol, lengthSymbol, timeSymbol);
return mus;
return ModelUnitSystem.createVCModelUnitSystem(volumeSubstanceSymbol, membraneSubstanceSymbol,
lumpedReactionSubstanceSymbol, volumeSymbol, areaSymbol, lengthSymbol, timeSymbol);
}
public static BioModel createBioModelWithSBMLUnitSystem(BioModel oldBioModel) throws ExpressionException, XmlParseException {
public static BioModel createBioModelWithSBMLUnitSystem(BioModel oldBioModel) throws ExpressionException, XmlParseException, ImageException, GeometryException {
ModelUnitSystem mus = createSbmlModelUnitSystem();
BioModel newBioModel = createBioModelWithNewUnitSystem(oldBioModel, mus);
return newBioModel;
return createBioModelWithNewUnitSystem(oldBioModel, mus);
}


public static BioModel createBioModelWithNewUnitSystem(BioModel oldBioModel, ModelUnitSystem newUnitSystem)
throws ExpressionException, XmlParseException {
throws ExpressionException, XmlParseException, ImageException, GeometryException {

oldBioModel.refreshDependencies();
Map<String, MathDescription> previousMathDescriptionMap = new LinkedHashMap<>();
Expand Down Expand Up @@ -75,12 +81,11 @@ public static BioModel createBioModelWithNewUnitSystem(BioModel oldBioModel, Mod

for (ReactionStep reactionStep : newBioModel.getModel().getReactionSteps()) {
SymbolTable oldSymbolTable = oldBioModel.getModel().getReactionStep(reactionStep.getName());
SymbolTable newSymbolTable = reactionStep;
for (Parameter p : reactionStep.getKinetics().getUnresolvedParameters()){
convertVarsWithUnitFactors(oldSymbolTable, newSymbolTable, p, dimensionless, KMOLE);
for (Parameter p : reactionStep.getKinetics().getUnresolvedParameters()){
convertVarsWithUnitFactors(oldSymbolTable, reactionStep, p, dimensionless, KMOLE);
}
for (Parameter p : reactionStep.getKinetics().getKineticsParameters()){
convertVarsWithUnitFactors(oldSymbolTable, newSymbolTable, p, dimensionless, KMOLE);
convertVarsWithUnitFactors(oldSymbolTable, reactionStep, p, dimensionless, KMOLE);
}

// We no longer have to deal with conversion factor expressions and try to be smart and simplify/flatten and rebind
Expand Down Expand Up @@ -173,7 +178,7 @@ public static BioModel createBioModelWithNewUnitSystem(BioModel oldBioModel, Mod
*/
// convert rate rules
RateRule[] rateRules = simContext.getRateRules();
if (rateRules != null && rateRules.length > 0) {
if (rateRules != null && rateRules.length > 0) {
for (RateRule rateRule : rateRules) {
RateRule oldRateRule = oldSimContext.getRateRule(rateRule.getName());
ScopedSymbolTable oldSymbolTable = oldRateRule.getSimulationContext();
Expand All @@ -198,6 +203,8 @@ public static BioModel createBioModelWithNewUnitSystem(BioModel oldBioModel, Mod
convertExprWithUnitFactors(oldSymbolTable, newSymbolTable, oldTargetUnit, newTargetUnit, assignmentRuleExpr, dimensionless, KMOLE);
}
}

convertGeometryWithUnitFactors(oldSimContext, simContext, dimensionless, KMOLE);
} // end for - simulationContext
newBioModel.refreshDependencies();
for (SimulationContext simContext : newBioModel.getSimulationContexts()){
Expand All @@ -215,28 +222,28 @@ public static BioModel createBioModelWithNewUnitSystem(BioModel oldBioModel, Mod
return newBioModel;
}

private static void convertVarsWithUnitFactors(SymbolTable oldSymbolTable, SymbolTable newSymbolTable, Parameter parameter,
private static void convertVarsWithUnitFactors(SymbolTable oldSymbolTable, SymbolTable newSymbolTable, Parameter newParameter,
VCUnitDefinition dimensionless, Model.ReservedSymbol KMOLE) throws ExpressionException {
// get old unit
VCUnitDefinition oldExprUnit = null;
Parameter oldParameter = (Parameter)oldSymbolTable.getEntry(parameter.getName());
Parameter oldParameter = (Parameter)oldSymbolTable.getEntry(newParameter.getName());
if (oldParameter == null){
System.err.println("parameter "+parameter.getName() + " was not found in the old symbol table");
System.err.println("parameter "+newParameter.getName() + " was not found in the old symbol table");
}else if (oldParameter.getUnitDefinition() == null){
System.err.println("parameter "+parameter.getName() + " has a null unit in old model, can't convert");
System.err.println("parameter "+newParameter.getName() + " has a null unit in old model, can't convert");
}else{
oldExprUnit = oldParameter.getUnitDefinition();
}

// get new unit
VCUnitDefinition newExprUnit = null;
if (parameter.getUnitDefinition() == null){
System.err.println("parameter "+parameter.getName() + " has a null unit in new model, can't convert");
if (newParameter.getUnitDefinition() == null){
System.err.println("parameter "+newParameter.getName() + " has a null unit in new model, can't convert");
}else{
newExprUnit = parameter.getUnitDefinition();
newExprUnit = newParameter.getUnitDefinition();
}

convertExprWithUnitFactors(oldSymbolTable, newSymbolTable, oldExprUnit, newExprUnit, parameter.getExpression(), dimensionless, KMOLE);
convertExprWithUnitFactors(oldSymbolTable, newSymbolTable, oldExprUnit, newExprUnit, newParameter.getExpression(), dimensionless, KMOLE);
}

/**
Expand Down Expand Up @@ -310,7 +317,6 @@ public static RationalExp getDimensionlessScaleFactorAsRationalExp(VCUnitDefinit
private static void convertExprWithUnitFactors(SymbolTable oldSymbolTable, SymbolTable newSymbolTable,
VCUnitDefinition oldExprUnit, VCUnitDefinition newExprUnit,
Expression expr, VCUnitDefinition dimensionless, Model.ReservedSymbol KMOLE) throws ExpressionException {

if (expr == null) {
return;
}
Expand Down Expand Up @@ -373,6 +379,76 @@ private static void convertExprWithUnitFactors(SymbolTable oldSymbolTable, Symbo
Expression origExp = new Expression(expr);
expr.substituteInPlace(origExp,flattened);
}

private static void convertGeometryWithUnitFactors(SimulationContext oldSimContext, SimulationContext newSimContext,
VCUnitDefinition dimensionless, Model.ReservedSymbol KMOLE)
throws ExpressionException, GeometryException, ImageException {
Geometry newGeometry = newSimContext.getGeometry();
VCUnitDefinition oldLengthUnit = oldSimContext.getModel().getUnitSystem().getLengthUnit();
VCUnitDefinition newLengthUnit = newSimContext.getModel().getUnitSystem().getLengthUnit();
VCUnitDefinition oldToNewConversionUnit = newLengthUnit.divideBy(oldLengthUnit);
Expression conversionFactor = getDimensionlessScaleFactor(oldToNewConversionUnit, dimensionless, KMOLE);

double lengthScaleFactor;
try {
lengthScaleFactor = ModelUnitConverter.round(conversionFactor.evaluateConstant(), 10);
} catch (Exception e){
throw new RuntimeException("Geometry scaling is unusual; failed to evaluate conversion factor", e);
}
if (lengthScaleFactor == 1) return;

// Scale Origin
Origin oldOrigin = newGeometry.getOrigin();
double correctedX = ModelUnitConverter.round(oldOrigin.getX() * lengthScaleFactor, 10);
if (correctedX == 0) correctedX = oldOrigin.getX() * lengthScaleFactor; // to prevent accidentally truncating to 0
double correctedY = ModelUnitConverter.round(oldOrigin.getY() * lengthScaleFactor, 10);
if (correctedY == 0) correctedY = oldOrigin.getY() * lengthScaleFactor; // to prevent accidentally truncating to 0
double correctedZ = ModelUnitConverter.round(oldOrigin.getZ() * lengthScaleFactor, 10);
if (correctedZ == 0) correctedZ = oldOrigin.getZ() * lengthScaleFactor; // to prevent accidentally truncating to 0
newGeometry.setOrigin(new Origin(correctedX, correctedY, correctedZ));

// Scale Extent
Extent oldExtent = newGeometry.getExtent();
correctedX = ModelUnitConverter.round(oldExtent.getX() * lengthScaleFactor, 10);
if (correctedX == 0) correctedX = oldExtent.getX() * lengthScaleFactor; // to prevent accidentally truncating to 0
correctedY = ModelUnitConverter.round(oldExtent.getY() * lengthScaleFactor, 10);
if (correctedY == 0) correctedY = oldExtent.getY() * lengthScaleFactor; // to prevent accidentally truncating to 0
correctedZ = ModelUnitConverter.round(oldExtent.getZ() * lengthScaleFactor, 10);
if (correctedZ == 0) correctedZ = oldExtent.getZ() * lengthScaleFactor; // to prevent accidentally truncating to 0
try {
newGeometry.setExtent(new Extent(correctedX, correctedY, correctedZ));
} catch (Exception e){
throw new RuntimeException("Problem with Extent: ", e);
}

// Scale analytical geometry
for (SubVolume subVol : newGeometry.getGeometrySpec().getSubVolumes()){
if (!(subVol instanceof AnalyticSubVolume anSubVar)) continue;
Expression newAnExpression = anSubVar.getExpression();
convertExprWithUnitFactors(oldSimContext, newSimContext, dimensionless, dimensionless, newAnExpression, dimensionless, KMOLE);
}

// Scale Constructed Solid Geometry
for (SubVolume subVol : newGeometry.getGeometrySpec().getSubVolumes()){
if (!(subVol instanceof CSGObject csg)) continue;
csg.rescaleInPlace(lengthScaleFactor);
}

GeometrySurfaceDescription gsd = newGeometry.getGeometrySurfaceDescription();
if (gsd != null){
newGeometry.getGeometrySpec().getSampledImage().setDirty();
gsd.updateAll();
}
newSimContext.refreshSpatialObjects();
}

private static double round(double value, int places) {
if (places < 0) throw new IllegalArgumentException();

BigDecimal bd = new BigDecimal(Double.toString(value));
bd = bd.setScale(places, RoundingMode.HALF_UP);
return bd.doubleValue();
}

public static void main(String[] args) {

Expand Down Expand Up @@ -429,10 +505,5 @@ public static void main(String[] args) {
System.out.println(a/b+"");
b = 1E-26;
System.out.println(a/b+"");





}
}
Expand Up @@ -10,7 +10,6 @@

package cbit.vcell.geometry;
import cbit.vcell.parser.*;
import org.vcell.util.Compare;
import org.vcell.util.Matchable;
import org.vcell.util.document.KeyValue;

Expand All @@ -24,7 +23,7 @@ public class AnalyticSubVolume extends SubVolume {
//
// the following fields are used for evaluating the inside/outside function "isInside()"
//
private double valueArray[] = new double[3];
private final double[] valueArray = new double[3];

/**
* MathDescription constructor comment.
Expand Down Expand Up @@ -70,19 +69,14 @@ public boolean compareEqual(Matchable obj) {
if (!compareEqual0(obj)){
return false;
}
if (!(obj instanceof AnalyticSubVolume)){
if (!(obj instanceof AnalyticSubVolume sv)){
return false;
}
AnalyticSubVolume sv = (AnalyticSubVolume)obj;

if ((exp==null && sv.exp!=null) || (exp!=null && sv.exp==null)){
if ((exp==null && sv.exp!=null) || (exp!=null && sv.exp==null)){
return false;
}
if (!ExpressionUtils.functionallyEquivalent(exp,sv.exp)){
return false;
}

return true;
return ExpressionUtils.functionallyEquivalent(exp, sv.exp);
}


Expand Down
Expand Up @@ -14,7 +14,6 @@

import cbit.vcell.render.Affine;

@SuppressWarnings("serial")
public class CSGHomogeneousTransformation extends CSGTransformation {

public CSGHomogeneousTransformation(String name, Affine forward, Affine inverse){
Expand All @@ -33,13 +32,8 @@ public boolean compareEqual(Matchable obj) {
if (!compareEqual0(obj)){
return false;
}
if (!(obj instanceof CSGHomogeneousTransformation)){
return false;
}
// CSGHomogeneousTransformation csght = (CSGHomogeneousTransformation)obj;

return true;
}
return obj instanceof CSGHomogeneousTransformation;
}

@Override
public CSGNode clone() {
Expand Down
15 changes: 5 additions & 10 deletions vcell-core/src/main/java/cbit/vcell/geometry/CSGNode.java
Expand Up @@ -17,7 +17,7 @@

import cbit.vcell.render.Vect3d;

@SuppressWarnings("serial")

public abstract class CSGNode implements Matchable, Serializable {
private String name;
protected CSGNode(String name) {
Expand All @@ -36,16 +36,11 @@ public final void setName(String name) {
}

protected boolean compareEqual0(Matchable obj) {
if (!(obj instanceof CSGNode)) {
return false;
}

CSGNode csgn = (CSGNode)obj;
if (!Compare.isEqual(name, csgn.name)){
if (!(obj instanceof CSGNode csgn)) {
return false;
}
return true;
}

return Compare.isEqual(name, csgn.name);
}

}
27 changes: 14 additions & 13 deletions vcell-core/src/main/java/cbit/vcell/geometry/CSGObject.java
Expand Up @@ -18,7 +18,6 @@
import cbit.vcell.parser.ExpressionException;
import cbit.vcell.render.Vect3d;

@SuppressWarnings("serial")
public class CSGObject extends SubVolume {

public static final String PROPERTY_NAME_ROOT = "root";
Expand All @@ -37,17 +36,12 @@ public boolean compareEqual(Matchable obj) {
if (!compareEqual0(obj)){
return false;
}
if (!(obj instanceof CSGObject)){
if (!(obj instanceof CSGObject csgo)){
return false;
}
CSGObject csgo = (CSGObject)obj;

if (!Compare.isEqualOrNull(root,csgo.root)){
return false;
}

return true;
}
return Compare.isEqualOrNull(root, csgo.root);
}

@Override
public boolean isInside(double x, double y, double z, GeometrySpec geometrySpec) throws GeometryException, ImageException, ExpressionException {
Expand Down Expand Up @@ -75,10 +69,7 @@ private CSGNode findCSGNodeByName(CSGNode node, String name) {
}
} else if (node instanceof CSGTransformation) {
CSGNode child = ((CSGTransformation) node).getChild();
CSGNode csgNode = findCSGNodeByName(child, name);
if (csgNode != null) {
return csgNode;
}
return findCSGNodeByName(child, name);
}
return null;
}
Expand Down Expand Up @@ -112,4 +103,14 @@ public void setRoot(CSGNode newValue) {
this.root = newValue;
firePropertyChange(PROPERTY_NAME_ROOT, oldValue, newValue);
}

public void rescaleInPlace(double lengthScaleFactor){
if (!(this.root instanceof CSGScale scaledNode)){
CSGNode temp = this.root;
this.root = new CSGScale("Unit Scale", new Vect3d(lengthScaleFactor, lengthScaleFactor, lengthScaleFactor));
((CSGScale)this.root).setChild(temp);
} else {
scaledNode.getScale().scale(lengthScaleFactor);
}
}
}

0 comments on commit 4c0c15b

Please sign in to comment.