use of cbit.vcell.mapping.StructureMapping.StructureMappingParameter in project vcell by virtualcell.
the class CurrentClampElectricalDevice method initializeParameters.
private void initializeParameters() throws ExpressionException {
ElectricalDevice.ElectricalDeviceParameter[] parameters = new ElectricalDevice.ElectricalDeviceParameter[3];
//
// set the transmembrane current (total current, if necessary derive it from the current density).
//
ElectricalDeviceParameter transMembraneCurrent = null;
ModelUnitSystem modelUnitSystem = mathMapping.getSimulationContext().getModel().getUnitSystem();
VCUnitDefinition currentUnit = modelUnitSystem.getCurrentUnit();
if (currentClampStimulus instanceof TotalCurrentClampStimulus) {
TotalCurrentClampStimulus stimulus = (TotalCurrentClampStimulus) currentClampStimulus;
LocalParameter currentParameter = stimulus.getCurrentParameter();
transMembraneCurrent = new ElectricalDeviceParameter(DefaultNames[ROLE_TransmembraneCurrent], new Expression(currentParameter.getExpression()), ROLE_TransmembraneCurrent, currentUnit);
} else if (currentClampStimulus instanceof CurrentDensityClampStimulus) {
CurrentDensityClampStimulus stimulus = (CurrentDensityClampStimulus) currentClampStimulus;
LocalParameter currentDensityParameter = stimulus.getCurrentDensityParameter();
//
// here we have to determine the expression for current (from current density).
//
Feature feature1 = currentClampStimulus.getElectrode().getFeature();
Feature feature2 = mathMapping.getSimulationContext().getGroundElectrode().getFeature();
Membrane membrane = null;
StructureTopology structTopology = mathMapping.getSimulationContext().getModel().getStructureTopology();
if (structTopology.getParentStructure(feature1) != null && structTopology.getOutsideFeature((Membrane) structTopology.getParentStructure(feature1)) == feature2) {
membrane = ((Membrane) structTopology.getParentStructure(feature1));
} else if (structTopology.getParentStructure(feature2) != null && structTopology.getOutsideFeature((Membrane) structTopology.getParentStructure(feature2)) == feature1) {
membrane = ((Membrane) structTopology.getParentStructure(feature2));
}
if (membrane == null) {
throw new RuntimeException("current clamp based on current density crosses multiple membranes, unable to " + "determine single membrane to convert current density into current in Application '" + mathMapping.getSimulationContext().getName() + "'.");
}
MembraneMapping membraneMapping = (MembraneMapping) mathMapping.getSimulationContext().getGeometryContext().getStructureMapping(membrane);
StructureMappingParameter sizeParameter = membraneMapping.getSizeParameter();
Expression area = null;
if (mathMapping.getSimulationContext().getGeometry().getDimension() == 0 && (sizeParameter.getExpression() == null || sizeParameter.getExpression().isZero())) {
area = membraneMapping.getNullSizeParameterValue();
} else {
area = new Expression(sizeParameter, mathMapping.getNameScope());
}
transMembraneCurrent = new ElectricalDeviceParameter(DefaultNames[ROLE_TransmembraneCurrent], Expression.mult(new Expression(currentDensityParameter.getExpression()), area), ROLE_TransmembraneCurrent, currentUnit);
} else {
throw new RuntimeException("unexpected current clamp stimulus type : " + currentClampStimulus.getClass().getName());
}
ElectricalDeviceParameter totalCurrent = new ElectricalDeviceParameter(DefaultNames[ROLE_TotalCurrent], new Expression(transMembraneCurrent, getNameScope()), ROLE_TotalCurrent, currentUnit);
ElectricalDeviceParameter voltage = new ElectricalDeviceParameter(DefaultNames[ROLE_Voltage], null, ROLE_Voltage, modelUnitSystem.getVoltageUnit());
parameters[0] = totalCurrent;
parameters[1] = transMembraneCurrent;
parameters[2] = voltage;
//
// add any user-defined parameters
//
LocalParameter[] stimulusParameters = currentClampStimulus.getLocalParameters();
for (int i = 0; stimulusParameters != null && i < stimulusParameters.length; i++) {
if (stimulusParameters[i].getRole() == ElectricalStimulus.ElectricalStimulusParameterType.UserDefined) {
ElectricalDeviceParameter newParam = new ElectricalDeviceParameter(stimulusParameters[i].getName(), new Expression(stimulusParameters[i].getExpression()), ROLE_UserDefined, stimulusParameters[i].getUnitDefinition());
parameters = (ElectricalDeviceParameter[]) BeanUtils.addElement(parameters, newParam);
}
}
setParameters(parameters);
}
use of cbit.vcell.mapping.StructureMapping.StructureMappingParameter in project vcell by virtualcell.
the class ElectricalCircuitGraph method getTotalMembraneCurrent.
/**
* Insert the method's description here.
* Creation date: (2/19/2002 12:56:13 PM)
* @return cbit.vcell.parser.Expression
* @param simContext cbit.vcell.mapping.SimulationContext
* @param membrane cbit.vcell.model.Membrane
*/
private static Expression getTotalMembraneCurrent(SimulationContext simContext, Membrane membrane, AbstractMathMapping mathMapping) throws ExpressionException {
MembraneMapping membraneMapping = (MembraneMapping) simContext.getGeometryContext().getStructureMapping(membrane);
if (!membraneMapping.getCalculateVoltage()) {
return new Expression(0.0);
}
//
// gather current terms
//
Expression currentExp = new Expression(0.0);
ReactionSpec[] reactionSpecs = simContext.getReactionContext().getReactionSpecs();
StructureMappingParameter sizeParameter = membraneMapping.getSizeParameter();
Expression area = null;
if (simContext.getGeometry().getDimension() == 0 && (sizeParameter.getExpression() == null || sizeParameter.getExpression().isZero())) {
System.out.println("size not set for membrane \"" + membrane.getName() + "\", refer to Structure Mapping in Application \"" + mathMapping.getSimulationContext().getName() + "\"");
area = membraneMapping.getNullSizeParameterValue();
} else {
area = new Expression(sizeParameter, mathMapping.getNameScope());
}
for (int i = 0; i < reactionSpecs.length; i++) {
//
if (reactionSpecs[i].isExcluded()) {
continue;
}
if (reactionSpecs[i].getReactionStep().getKinetics() instanceof DistributedKinetics) {
ReactionStep rs = reactionSpecs[i].getReactionStep();
DistributedKinetics distributedKinetics = (DistributedKinetics) rs.getKinetics();
if (rs.getStructure() == membrane) {
if (!distributedKinetics.getCurrentDensityParameter().getExpression().isZero()) {
//
// change sign convension from inward current to outward current (which is consistent with voltage convension)
//
currentExp = Expression.add(currentExp, Expression.negate(Expression.mult(new Expression(distributedKinetics.getCurrentDensityParameter(), mathMapping.getNameScope()), area)));
}
}
} else {
ReactionStep rs = reactionSpecs[i].getReactionStep();
LumpedKinetics lumpedKinetics = (LumpedKinetics) rs.getKinetics();
if (rs.getStructure() == membrane) {
if (!lumpedKinetics.getLumpedCurrentParameter().getExpression().isZero()) {
//
if (!(membraneMapping.getGeometryClass() instanceof CompartmentSubVolume)) {
throw new RuntimeException("math generation for total currents within spatial electrophysiology not yet implemented");
}
Expression lumpedCurrentSymbolExp = new Expression(lumpedKinetics.getLumpedCurrentParameter(), mathMapping.getNameScope());
currentExp = Expression.add(currentExp, Expression.negate(lumpedCurrentSymbolExp));
}
}
}
}
return currentExp.flatten();
}
use of cbit.vcell.mapping.StructureMapping.StructureMappingParameter in project vcell by virtualcell.
the class CurrentClampElectricalDevice method initializeParameters.
private void initializeParameters() throws ExpressionException {
ElectricalDevice.ElectricalDeviceParameter[] parameters = new ElectricalDevice.ElectricalDeviceParameter[3];
//
// set the transmembrane current (total current, if necessary derive it from the current density).
//
ElectricalDeviceParameter transMembraneCurrent = null;
ModelUnitSystem modelUnitSystem = mathMapping_4_8.getSimulationContext().getModel().getUnitSystem();
VCUnitDefinition currentUnit = modelUnitSystem.getCurrentUnit();
if (currentClampStimulus instanceof TotalCurrentClampStimulus) {
TotalCurrentClampStimulus stimulus = (TotalCurrentClampStimulus) currentClampStimulus;
LocalParameter currentParameter = stimulus.getCurrentParameter();
transMembraneCurrent = new ElectricalDeviceParameter(DefaultNames[ROLE_TransmembraneCurrent], new Expression(currentParameter.getExpression()), ROLE_TransmembraneCurrent, currentUnit);
} else if (currentClampStimulus instanceof CurrentDensityClampStimulus) {
CurrentDensityClampStimulus stimulus = (CurrentDensityClampStimulus) currentClampStimulus;
LocalParameter currentDensityParameter = stimulus.getCurrentDensityParameter();
//
// here we have to determine the expression for current (from current density).
//
Feature feature1 = currentClampStimulus.getElectrode().getFeature();
Feature feature2 = mathMapping_4_8.getSimulationContext().getGroundElectrode().getFeature();
Membrane membrane = null;
StructureTopology structTopology = mathMapping_4_8.getSimulationContext().getModel().getStructureTopology();
if (structTopology.getParentStructure(feature1) != null && structTopology.getOutsideFeature((Membrane) structTopology.getParentStructure(feature1)) == feature2) {
membrane = ((Membrane) structTopology.getParentStructure(feature1));
} else if (structTopology.getParentStructure(feature2) != null && structTopology.getOutsideFeature((Membrane) structTopology.getParentStructure(feature2)) == feature1) {
membrane = ((Membrane) structTopology.getParentStructure(feature2));
}
if (membrane == null) {
throw new RuntimeException("current clamp based on current density crosses multiple membranes, unable to " + "determine single membrane to convert current density into current in Application '" + mathMapping_4_8.getSimulationContext().getName() + "'.");
}
MembraneMapping membraneMapping = (MembraneMapping) mathMapping_4_8.getSimulationContext().getGeometryContext().getStructureMapping(membrane);
StructureMappingParameter sizeParameter = membraneMapping.getSizeParameter();
Expression area = null;
if (mathMapping_4_8.getSimulationContext().getGeometry().getDimension() == 0 && (sizeParameter.getExpression() == null || sizeParameter.getExpression().isZero())) {
area = membraneMapping.getNullSizeParameterValue();
} else {
area = new Expression(sizeParameter, mathMapping_4_8.getNameScope());
}
transMembraneCurrent = new ElectricalDeviceParameter(DefaultNames[ROLE_TransmembraneCurrent], Expression.mult(new Expression(currentDensityParameter.getExpression()), area), ROLE_TransmembraneCurrent, currentUnit);
} else {
throw new RuntimeException("unexpected current clamp stimulus type : " + currentClampStimulus.getClass().getName());
}
ElectricalDeviceParameter totalCurrent = new ElectricalDeviceParameter(DefaultNames[ROLE_TotalCurrent], new Expression(transMembraneCurrent, getNameScope()), ROLE_TotalCurrent, currentUnit);
ElectricalDeviceParameter voltage = new ElectricalDeviceParameter(DefaultNames[ROLE_Voltage], null, ROLE_Voltage, modelUnitSystem.getVoltageUnit());
parameters[0] = totalCurrent;
parameters[1] = transMembraneCurrent;
parameters[2] = voltage;
//
// add any user-defined parameters
//
LocalParameter[] stimulusParameters = currentClampStimulus.getLocalParameters();
for (int i = 0; stimulusParameters != null && i < stimulusParameters.length; i++) {
if (stimulusParameters[i].getRole() == ElectricalStimulus.ElectricalStimulusParameterType.UserDefined) {
ElectricalDeviceParameter newParam = new ElectricalDeviceParameter(stimulusParameters[i].getName(), new Expression(stimulusParameters[i].getExpression()), ROLE_UserDefined, stimulusParameters[i].getUnitDefinition());
parameters = (ElectricalDeviceParameter[]) BeanUtils.addElement(parameters, newParam);
}
}
setParameters(parameters);
}
use of cbit.vcell.mapping.StructureMapping.StructureMappingParameter in project vcell by virtualcell.
the class GeometryContext method setDefaultUnitSizes.
private void setDefaultUnitSizes() throws PropertyVetoException {
if (getGeometry().getDimension() > 0) {
for (StructureMapping sm : fieldStructureMappings) {
StructureMapping[] sms = getStructureMappings(sm.getGeometryClass());
StructureMappingParameter unitSizeParameter = sm.getUnitSizeParameter();
if (unitSizeParameter == null) {
continue;
}
Expression exp = unitSizeParameter.getExpression();
if (sm instanceof MembraneMapping) {
// Membrane mapped to surface or subdomain, default to 1.0
if (exp == null) {
try {
unitSizeParameter.setExpression(new Expression(1.0));
} catch (ExpressionBindingException e) {
e.printStackTrace();
}
}
} else if (sm instanceof FeatureMapping) {
// Feature mapped to subdomain
if (sm.getGeometryClass() instanceof SubVolume) {
if (sms != null && sms.length == 1) {
try {
unitSizeParameter.setExpression(new Expression(1.0));
} catch (ExpressionBindingException e) {
e.printStackTrace();
}
}
} else {
if (exp == null) {
try {
unitSizeParameter.setExpression(new Expression(1.0));
} catch (ExpressionBindingException e) {
e.printStackTrace();
}
}
}
}
}
}
}
use of cbit.vcell.mapping.StructureMapping.StructureMappingParameter in project vcell by virtualcell.
the class SBMLExporter method addGeometry.
private void addGeometry() throws SbmlException {
SpatialModelPlugin mplugin = (SpatialModelPlugin) sbmlModel.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
// Creates a geometry object via SpatialModelPlugin object.
org.sbml.jsbml.ext.spatial.Geometry sbmlGeometry = mplugin.createGeometry();
sbmlGeometry.setCoordinateSystem(GeometryKind.cartesian);
sbmlGeometry.setSpatialId("vcell");
Geometry vcGeometry = getSelectedSimContext().getGeometry();
Model vcModel = getSelectedSimContext().getModel();
//
// list of CoordinateComponents : 1 if geometry is 1-d, 2 if geometry is 2-d, 3 if geometry is 3-d
//
int dimension = vcGeometry.getDimension();
Extent vcExtent = vcGeometry.getExtent();
Origin vcOrigin = vcGeometry.getOrigin();
// add x coordinate component
CoordinateComponent xComp = sbmlGeometry.createCoordinateComponent();
xComp.setSpatialId(vcModel.getX().getName());
xComp.setType(CoordinateKind.cartesianX);
final UnitDefinition sbmlUnitDef_length = getOrCreateSBMLUnit(vcModel.getUnitSystem().getLengthUnit());
xComp.setUnits(sbmlUnitDef_length);
Boundary minX = new Boundary();
xComp.setBoundaryMinimum(minX);
minX.setSpatialId("Xmin");
minX.setValue(vcOrigin.getX());
Boundary maxX = new Boundary();
xComp.setBoundaryMaximum(maxX);
maxX.setSpatialId("Xmax");
maxX.setValue(vcOrigin.getX() + (vcExtent.getX()));
org.sbml.jsbml.Parameter pX = sbmlModel.createParameter();
pX.setId(vcModel.getX().getName());
pX.setValue(0.0);
pX.setConstant(false);
pX.setUnits(sbmlUnitDef_length);
SpatialParameterPlugin spPluginPx = (SpatialParameterPlugin) pX.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
SpatialSymbolReference spSymRefPx = new SpatialSymbolReference();
spPluginPx.setParamType(spSymRefPx);
spSymRefPx.setSpatialRef(xComp.getSpatialId());
// add y coordinate component
if (dimension == 2 || dimension == 3) {
CoordinateComponent yComp = sbmlGeometry.createCoordinateComponent();
yComp.setSpatialId(vcModel.getY().getName());
yComp.setType(CoordinateKind.cartesianY);
yComp.setUnits(sbmlUnitDef_length);
Boundary minY = new Boundary();
yComp.setBoundaryMinimum(minY);
minY.setSpatialId("Ymin");
minY.setValue(vcOrigin.getY());
Boundary maxY = new Boundary();
yComp.setBoundaryMaximum(maxY);
maxY.setSpatialId("Ymax");
maxY.setValue(vcOrigin.getY() + (vcExtent.getY()));
org.sbml.jsbml.Parameter pY = sbmlModel.createParameter();
pY.setId(vcModel.getY().getName());
pY.setValue(0.0);
pY.setConstant(false);
pY.setUnits(sbmlUnitDef_length);
SpatialParameterPlugin spPluginPy = (SpatialParameterPlugin) pY.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
SpatialSymbolReference spSymRefPy = new SpatialSymbolReference();
spPluginPy.setParamType(spSymRefPy);
spSymRefPy.setSpatialRef(yComp.getSpatialId());
}
// add z coordinate component
if (dimension == 3) {
CoordinateComponent zComp = sbmlGeometry.createCoordinateComponent();
zComp.setSpatialId(vcModel.getZ().getName());
zComp.setType(CoordinateKind.cartesianZ);
zComp.setUnits(sbmlUnitDef_length);
Boundary minZ = new Boundary();
zComp.setBoundaryMinimum(minZ);
minZ.setSpatialId("Zmin");
minZ.setValue(vcOrigin.getZ());
Boundary maxZ = new Boundary();
zComp.setBoundaryMaximum(maxZ);
maxZ.setSpatialId("Zmax");
maxZ.setValue(vcOrigin.getZ() + (vcExtent.getZ()));
org.sbml.jsbml.Parameter pZ = sbmlModel.createParameter();
pZ.setId(vcModel.getZ().getName());
pZ.setValue(0.0);
pZ.setConstant(false);
pZ.setUnits(sbmlUnitDef_length);
SpatialParameterPlugin spPluginPz = (SpatialParameterPlugin) pZ.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
SpatialSymbolReference spSymRefPz = new SpatialSymbolReference();
spPluginPz.setParamType(spSymRefPz);
spSymRefPz.setSpatialRef(zComp.getSpatialId());
}
//
// list of compartmentMappings : VC structureMappings
//
GeometryContext vcGeoContext = getSelectedSimContext().getGeometryContext();
StructureMapping[] vcStrucMappings = vcGeoContext.getStructureMappings();
for (int i = 0; i < vcStrucMappings.length; i++) {
StructureMapping vcStructMapping = vcStrucMappings[i];
String structName = vcStructMapping.getStructure().getName();
Compartment comp = sbmlModel.getCompartment(TokenMangler.mangleToSName(structName));
SpatialCompartmentPlugin cplugin = (SpatialCompartmentPlugin) comp.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
GeometryClass gc = vcStructMapping.getGeometryClass();
if (!goodPointer(gc, GeometryClass.class, structName)) {
continue;
}
CompartmentMapping compMapping = new CompartmentMapping();
cplugin.setCompartmentMapping(compMapping);
String geomClassName = gc.getName();
String id = TokenMangler.mangleToSName(geomClassName + structName);
compMapping.setSpatialId(id);
compMapping.setDomainType(TokenMangler.mangleToSName(DOMAIN_TYPE_PREFIX + geomClassName));
try {
StructureMappingParameter usp = vcStructMapping.getUnitSizeParameter();
Expression e = usp.getExpression();
if (goodPointer(e, Expression.class, id)) {
compMapping.setUnitSize(e.evaluateConstant());
}
} catch (ExpressionException e) {
e.printStackTrace(System.out);
throw new RuntimeException("Unable to create compartment mapping for structureMapping '" + compMapping.getId() + "' : " + e.getMessage());
}
}
//
// list of domain types : subvolumes and surface classes from VC
//
boolean bAnyAnalyticSubvolumes = false;
boolean bAnyImageSubvolumes = false;
boolean bAnyCSGSubvolumes = false;
GeometryClass[] vcGeomClasses = vcGeometry.getGeometryClasses();
int numSubVols = 0;
for (int i = 0; i < vcGeomClasses.length; i++) {
DomainType domainType = sbmlGeometry.createDomainType();
domainType.setSpatialId(DOMAIN_TYPE_PREFIX + vcGeomClasses[i].getName());
if (vcGeomClasses[i] instanceof SubVolume) {
if (((SubVolume) vcGeomClasses[i]) instanceof AnalyticSubVolume) {
bAnyAnalyticSubvolumes = true;
} else if (((SubVolume) vcGeomClasses[i]) instanceof ImageSubVolume) {
bAnyImageSubvolumes = true;
} else if (((SubVolume) vcGeomClasses[i]) instanceof CSGObject) {
bAnyCSGSubvolumes = true;
}
domainType.setSpatialDimensions(3);
numSubVols++;
} else if (vcGeomClasses[i] instanceof SurfaceClass) {
domainType.setSpatialDimensions(2);
}
}
//
// list of domains, adjacent domains : from VC geometricRegions
//
GeometrySurfaceDescription vcGSD = vcGeometry.getGeometrySurfaceDescription();
if (vcGSD.getRegionImage() == null) {
try {
vcGSD.updateAll();
} catch (Exception e) {
e.printStackTrace(System.out);
throw new RuntimeException("Unable to generate region images for geometry");
}
}
GeometricRegion[] vcGeometricRegions = vcGSD.getGeometricRegions();
ISize sampleSize = vcGSD.getVolumeSampleSize();
int numX = sampleSize.getX();
int numY = sampleSize.getY();
int numZ = sampleSize.getZ();
double ox = vcOrigin.getX();
double oy = vcOrigin.getY();
double oz = vcOrigin.getZ();
RegionInfo[] regionInfos = vcGSD.getRegionImage().getRegionInfos();
for (int i = 0; i < vcGeometricRegions.length; i++) {
// domains
Domain domain = sbmlGeometry.createDomain();
domain.setSpatialId(vcGeometricRegions[i].getName());
if (vcGeometricRegions[i] instanceof VolumeGeometricRegion) {
domain.setDomainType(DOMAIN_TYPE_PREFIX + ((VolumeGeometricRegion) vcGeometricRegions[i]).getSubVolume().getName());
//
// get a list of interior points ... should probably use the distance map to find a point
// furthest inside (or several points associated with the morphological skeleton).
//
InteriorPoint interiorPt = domain.createInteriorPoint();
int regionID = ((VolumeGeometricRegion) vcGeometricRegions[i]).getRegionID();
boolean bFound = false;
int regInfoIndx = 0;
for (int j = 0; j < regionInfos.length; j++) {
regInfoIndx = j;
if (regionInfos[j].getRegionIndex() == regionID) {
int volIndx = 0;
for (int z = 0; z < numZ && !bFound; z++) {
for (int y = 0; y < numY && !bFound; y++) {
for (int x = 0; x < numX && !bFound; x++) {
if (regionInfos[j].isIndexInRegion(volIndx)) {
bFound = true;
double unit_z = (numZ > 1) ? ((double) z) / (numZ - 1) : 0.5;
double coordZ = oz + vcExtent.getZ() * unit_z;
double unit_y = (numY > 1) ? ((double) y) / (numY - 1) : 0.5;
double coordY = oy + vcExtent.getY() * unit_y;
double unit_x = (numX > 1) ? ((double) x) / (numX - 1) : 0.5;
double coordX = ox + vcExtent.getX() * unit_x;
interiorPt.setCoord1(coordX);
interiorPt.setCoord2(coordY);
interiorPt.setCoord3(coordZ);
}
volIndx++;
}
// end - for x
}
// end - for y
}
// end - for z
}
// end if
}
// end for regionInfos
if (!bFound) {
throw new RuntimeException("Unable to find interior point for region '" + regionInfos[regInfoIndx].toString());
}
} else if (vcGeometricRegions[i] instanceof SurfaceGeometricRegion) {
SurfaceGeometricRegion vcSurfaceGeomReg = (SurfaceGeometricRegion) vcGeometricRegions[i];
GeometricRegion geomRegion0 = vcSurfaceGeomReg.getAdjacentGeometricRegions()[0];
GeometricRegion geomRegion1 = vcSurfaceGeomReg.getAdjacentGeometricRegions()[1];
SurfaceClass surfaceClass = vcGSD.getSurfaceClass(((VolumeGeometricRegion) geomRegion0).getSubVolume(), ((VolumeGeometricRegion) geomRegion1).getSubVolume());
domain.setDomainType(DOMAIN_TYPE_PREFIX + surfaceClass.getName());
// adjacent domains : 2 adjacent domain objects for each surfaceClass in VC.
// adjacent domain 1
GeometricRegion adjGeomRegion0 = vcSurfaceGeomReg.getAdjacentGeometricRegions()[0];
GeometricRegion adjGeomRegion1 = vcSurfaceGeomReg.getAdjacentGeometricRegions()[1];
AdjacentDomains adjDomain = new AdjacentDomains();
adjDomain.setSpatialId(TokenMangler.mangleToSName(vcSurfaceGeomReg.getName() + "_" + adjGeomRegion0.getName()));
adjDomain.setDomain1(vcSurfaceGeomReg.getName());
adjDomain.setDomain2(adjGeomRegion0.getName());
sbmlGeometry.addAdjacentDomain(adjDomain);
// adj domain 2
adjDomain = new AdjacentDomains();
adjDomain.setSpatialId(TokenMangler.mangleToSName(vcSurfaceGeomReg.getName() + "_" + adjGeomRegion1.getName()));
adjDomain.setDomain1(vcSurfaceGeomReg.getName());
adjDomain.setDomain2(adjGeomRegion1.getName());
sbmlGeometry.addAdjacentDomain(adjDomain);
}
}
//
if (bAnyAnalyticSubvolumes && !bAnyImageSubvolumes && !bAnyCSGSubvolumes) {
AnalyticGeometry sbmlAnalyticGeomDefinition = sbmlGeometry.createAnalyticGeometry();
sbmlAnalyticGeomDefinition.setSpatialId(TokenMangler.mangleToSName("Analytic_" + vcGeometry.getName()));
sbmlAnalyticGeomDefinition.setIsActive(true);
for (int i = 0; i < vcGeomClasses.length; i++) {
if (vcGeomClasses[i] instanceof AnalyticSubVolume) {
AnalyticVolume analyticVol = sbmlAnalyticGeomDefinition.createAnalyticVolume();
analyticVol.setSpatialId(vcGeomClasses[i].getName());
analyticVol.setDomainType(DOMAIN_TYPE_PREFIX + vcGeomClasses[i].getName());
analyticVol.setFunctionType(FunctionKind.layered);
analyticVol.setOrdinal(numSubVols - (i + 1));
Expression expr = ((AnalyticSubVolume) vcGeomClasses[i]).getExpression();
try {
String mathMLStr = ExpressionMathMLPrinter.getMathML(expr, true, MathType.BOOLEAN);
ASTNode mathMLNode = ASTNode.readMathMLFromString(mathMLStr);
analyticVol.setMath(mathMLNode);
} catch (Exception e) {
e.printStackTrace(System.out);
throw new RuntimeException("Error converting VC subvolume expression to mathML" + e.getMessage());
}
}
}
}
//
if (!bAnyAnalyticSubvolumes && !bAnyImageSubvolumes && bAnyCSGSubvolumes) {
CSGeometry sbmlCSGeomDefinition = new CSGeometry();
sbmlGeometry.addGeometryDefinition(sbmlCSGeomDefinition);
sbmlCSGeomDefinition.setSpatialId(TokenMangler.mangleToSName("CSG_" + vcGeometry.getName()));
for (int i = 0; i < vcGeomClasses.length; i++) {
if (vcGeomClasses[i] instanceof CSGObject) {
CSGObject vcellCSGObject = (CSGObject) vcGeomClasses[i];
org.sbml.jsbml.ext.spatial.CSGObject sbmlCSGObject = new org.sbml.jsbml.ext.spatial.CSGObject();
sbmlCSGeomDefinition.addCSGObject(sbmlCSGObject);
sbmlCSGObject.setSpatialId(vcellCSGObject.getName());
sbmlCSGObject.setDomainType(DOMAIN_TYPE_PREFIX + vcellCSGObject.getName());
// the ordinal should the the least for the default/background subVolume
sbmlCSGObject.setOrdinal(numSubVols - (i + 1));
org.sbml.jsbml.ext.spatial.CSGNode sbmlcsgNode = getSBMLCSGNode(vcellCSGObject.getRoot());
sbmlCSGObject.setCSGNode(sbmlcsgNode);
}
}
}
//
// add "Segmented" and "DistanceMap" SampledField Geometries
//
final boolean bVCGeometryIsImage = bAnyImageSubvolumes && !bAnyAnalyticSubvolumes && !bAnyCSGSubvolumes;
// 55if (bAnyAnalyticSubvolumes || bAnyImageSubvolumes || bAnyCSGSubvolumes){
if (bVCGeometryIsImage) {
//
// add "Segmented" SampledFieldGeometry
//
SampledFieldGeometry segmentedImageSampledFieldGeometry = sbmlGeometry.createSampledFieldGeometry();
segmentedImageSampledFieldGeometry.setSpatialId(TokenMangler.mangleToSName("SegmentedImage_" + vcGeometry.getName()));
segmentedImageSampledFieldGeometry.setIsActive(true);
// 55boolean bVCGeometryIsImage = bAnyImageSubvolumes && !bAnyAnalyticSubvolumes && !bAnyCSGSubvolumes;
Geometry vcImageGeometry = null;
{
if (bVCGeometryIsImage) {
// make a resampled image;
if (dimension == 3) {
try {
ISize imageSize = vcGeometry.getGeometrySpec().getDefaultSampledImageSize();
vcGeometry.precomputeAll(new GeometryThumbnailImageFactoryAWT());
vcImageGeometry = RayCaster.resampleGeometry(new GeometryThumbnailImageFactoryAWT(), vcGeometry, imageSize);
} catch (Throwable e) {
e.printStackTrace(System.out);
throw new RuntimeException("Unable to convert the original analytic or constructed solid geometry to image-based geometry : " + e.getMessage());
}
} else {
try {
vcGeometry.precomputeAll(new GeometryThumbnailImageFactoryAWT(), true, false);
GeometrySpec origGeometrySpec = vcGeometry.getGeometrySpec();
VCImage newVCImage = origGeometrySpec.getSampledImage().getCurrentValue();
//
// construct the new geometry with the sampled VCImage.
//
vcImageGeometry = new Geometry(vcGeometry.getName() + "_asImage", newVCImage);
vcImageGeometry.getGeometrySpec().setExtent(vcGeometry.getExtent());
vcImageGeometry.getGeometrySpec().setOrigin(vcGeometry.getOrigin());
vcImageGeometry.setDescription(vcGeometry.getDescription());
vcImageGeometry.getGeometrySurfaceDescription().setFilterCutoffFrequency(vcGeometry.getGeometrySurfaceDescription().getFilterCutoffFrequency());
vcImageGeometry.precomputeAll(new GeometryThumbnailImageFactoryAWT(), true, true);
} catch (Exception e) {
e.printStackTrace(System.out);
throw new RuntimeException("Unable to convert the original analytic or constructed solid geometry to image-based geometry : " + e.getMessage());
}
}
GeometryClass[] vcImageGeomClasses = vcImageGeometry.getGeometryClasses();
for (int j = 0; j < vcImageGeomClasses.length; j++) {
if (vcImageGeomClasses[j] instanceof ImageSubVolume) {
SampledVolume sampledVol = segmentedImageSampledFieldGeometry.createSampledVolume();
sampledVol.setSpatialId(vcGeomClasses[j].getName());
sampledVol.setDomainType(DOMAIN_TYPE_PREFIX + vcGeomClasses[j].getName());
sampledVol.setSampledValue(((ImageSubVolume) vcImageGeomClasses[j]).getPixelValue());
}
}
// add sampledField to sampledFieldGeometry
SampledField segmentedImageSampledField = sbmlGeometry.createSampledField();
VCImage vcImage = vcImageGeometry.getGeometrySpec().getImage();
segmentedImageSampledField.setSpatialId("SegmentedImageSampledField");
segmentedImageSampledField.setNumSamples1(vcImage.getNumX());
segmentedImageSampledField.setNumSamples2(vcImage.getNumY());
segmentedImageSampledField.setNumSamples3(vcImage.getNumZ());
segmentedImageSampledField.setInterpolationType(InterpolationKind.nearestneighbor);
segmentedImageSampledField.setCompression(CompressionKind.uncompressed);
segmentedImageSampledField.setDataType(DataKind.UINT8);
segmentedImageSampledFieldGeometry.setSampledField(segmentedImageSampledField.getId());
try {
byte[] vcImagePixelsBytes = vcImage.getPixels();
// imageData.setCompression("");
StringBuffer sb = new StringBuffer();
for (int i = 0; i < vcImagePixelsBytes.length; i++) {
int uint8_sample = ((int) vcImagePixelsBytes[i]) & 0xff;
sb.append(uint8_sample + " ");
}
segmentedImageSampledField.setSamplesLength(vcImage.getNumXYZ());
segmentedImageSampledField.setSamples(sb.toString().trim());
} catch (ImageException e) {
e.printStackTrace(System.out);
throw new RuntimeException("Unable to export image from VCell to SBML : " + e.getMessage());
}
}
}
/*
//
// add "DistanceMap" SampledFieldGeometry if there are exactly two subvolumes (else need more fields) and geometry is 3d.
//
if (numSubVols==2 && dimension == 3){
SignedDistanceMap[] distanceMaps = null;
try {
distanceMaps = DistanceMapGenerator.computeDistanceMaps(vcImageGeometry, vcImageGeometry.getGeometrySpec().getImage(), false, false);
} catch (ImageException e) {
e.printStackTrace(System.out);
System.err.println("Unable to export distance map sampled field from VCell to SBML : " + e.getMessage());
// throw new RuntimeException("Unable to export distance map sampled field from VCell to SBML : " + e.getMessage());
// don't want to throw an exception and stop export because distance map geometry couldn't be exported.
// just 'return' from method (since this is the last thing that is being done in this method).
return;
}
//
// the two distanceMaps should be redundant (one is negation of the other) ... so choose first one for field.
//
double[] signedDistances = distanceMaps[0].getSignedDistances();
SampledFieldGeometry distanceMapSampledFieldGeometry = sbmlGeometry.createSampledFieldGeometry();
distanceMapSampledFieldGeometry.setSpatialId(TokenMangler.mangleToSName("DistanceMap_"+vcGeometry.getName()));
SampledField distanceMapSampledField = distanceMapSampledFieldGeometry.createSampledField();
distanceMapSampledField.setSpatialId("DistanceMapSampledField");
distanceMapSampledField.setNumSamples1(distanceMaps[0].getSamplesX().length);
distanceMapSampledField.setNumSamples2(distanceMaps[0].getSamplesY().length);
distanceMapSampledField.setNumSamples3(distanceMaps[0].getSamplesZ().length);
distanceMapSampledField.setDataType("real");
System.err.println("do we need distanceMapSampleField.setDataType()?");
distanceMapSampledField.setInterpolationType("linear");
ImageData distanceMapImageData = distanceMapSampledField.createImageData();
distanceMapImageData.setDataType("int16");
System.err.println("should be:\n distanceMapImageData.setDataType(\"float32\")");
// distanceMapImageData.setCompression("");
double maxAbsValue = 0;
for (int i = 0; i < signedDistances.length; i++) {
maxAbsValue = Math.max(maxAbsValue,Math.abs(signedDistances[i]));
}
if (maxAbsValue==0.0){
throw new RuntimeException("computed distance map all zeros");
}
double scale = (Short.MAX_VALUE-1)/maxAbsValue;
int[] scaledIntegerDistanceMap = new int[signedDistances.length];
for (int i = 0; i < signedDistances.length; i++) {
scaledIntegerDistanceMap[i] = (int)(scale * signedDistances[i]);
}
distanceMapImageData.setSamples(scaledIntegerDistanceMap, signedDistances.length);
System.err.println("should be:\n distanceMapImageData.setSamples((float[])signedDistances,signedDistances.length)");
SampledVolume sampledVol = distanceMapSampledFieldGeometry.createSampledVolume();
sampledVol.setSpatialId(distanceMaps[0].getInsideSubvolumeName());
sampledVol.setDomainType(DOMAIN_TYPE_PREFIX+distanceMaps[0].getInsideSubvolumeName());
sampledVol.setSampledValue(255);
sampledVol = distanceMapSampledFieldGeometry.createSampledVolume();
sampledVol.setSpatialId(distanceMaps[1].getInsideSubvolumeName());
sampledVol.setDomainType(DOMAIN_TYPE_PREFIX+distanceMaps[1].getInsideSubvolumeName());
sampledVol.setSampledValue(1);
}
*/
}
//
// add "SurfaceMesh" ParametricGeometry
//
// if (bAnyAnalyticSubvolumes || bAnyImageSubvolumes || bAnyCSGSubvolumes){
// ParametricGeometry sbmlParametricGeomDefinition = sbmlGeometry.createParametricGeometry();
// sbmlParametricGeomDefinition.setSpatialId(TokenMangler.mangleToSName("SurfaceMesh_"+vcGeometry.getName()));
// xxxx
// }
}
Aggregations