Examples with SBMLImporter org.vcell.sbml.vcell.SBMLImporter used on opensource projects

Example 6 with SBMLImporter

use of org.vcell.sbml.vcell.SBMLImporter in project vcell by virtualcell.

the class VCellSBMLSolver method solveVCell.

public File solveVCell(String filePrefix, File outDir, String sbmlFileName, SimSpec testSpec) throws IOException, SolverException, SbmlException {
    try {
        cbit.util.xml.VCLogger logger = new LocalLogger();
        // 
        // Instantiate an SBMLImporter to get the speciesUnitsHash - to compute the conversion factor from VC->SB species units.
        // and import SBML  (sbml->bioModel)
        org.vcell.sbml.vcell.SBMLImporter sbmlImporter = new org.vcell.sbml.vcell.SBMLImporter(sbmlFileName, logger, false);
        BioModel bioModel = sbmlImporter.getBioModel();
        if (bRoundTrip) {
            // Round trip the bioModel (bioModel->sbml->bioModel).
            // export bioModel as sbml and save
            String vcml_sbml = cbit.vcell.xml.XmlHelper.exportSBML(bioModel, 2, 1, 0, false, bioModel.getSimulationContext(0), null);
            // re-import bioModel from exported sbml
            XMLSource vcml_sbml_Src = new XMLSource(vcml_sbml);
            BioModel newBioModel = (BioModel) XmlHelper.importSBML(logger, vcml_sbml_Src, false);
            // have rest of code use the round-tripped biomodel
            bioModel = newBioModel;
        }
        // 
        // select only Application, generate math, and create a single Simulation.
        // 
        SimulationContext simContext = bioModel.getSimulationContext(0);
        MathMapping mathMapping = simContext.createNewMathMapping();
        MathDescription mathDesc = mathMapping.getMathDescription();
        simContext.setMathDescription(mathDesc);
        SimulationVersion simVersion = new SimulationVersion(new KeyValue("100"), "unnamed", null, null, null, null, null, null, null, null);
        Simulation sim = new Simulation(simVersion, mathDesc);
        sim.setName("unnamed");
        // if time factor from SBML is not 1 (i.e., it is not in secs but in minutes or hours), convert endTime to min/hr as : endTime*timeFactor
        // double endTime = testSpec.getEndTime()*timeFactor;
        double endTime = testSpec.getEndTime();
        sim.getSolverTaskDescription().setTimeBounds(new TimeBounds(0, endTime));
        TimeStep timeStep = new TimeStep();
        sim.getSolverTaskDescription().setTimeStep(new TimeStep(timeStep.getMinimumTimeStep(), timeStep.getDefaultTimeStep(), endTime / 10000));
        sim.getSolverTaskDescription().setOutputTimeSpec(new UniformOutputTimeSpec((endTime - 0) / testSpec.getNumTimeSteps()));
        sim.getSolverTaskDescription().setErrorTolerance(new ErrorTolerance(testSpec.getAbsTolerance(), testSpec.getRelTolerance()));
        // sim.getSolverTaskDescription().setErrorTolerance(new ErrorTolerance(1e-10, 1e-12));
        // Generate .idaInput string
        File idaInputFile = new File(outDir, filePrefix + SimDataConstants.IDAINPUT_DATA_EXTENSION);
        PrintWriter idaPW = new java.io.PrintWriter(idaInputFile);
        SimulationJob simJob = new SimulationJob(sim, 0, null);
        SimulationTask simTask = new SimulationTask(simJob, 0);
        IDAFileWriter idaFileWriter = new IDAFileWriter(idaPW, simTask);
        idaFileWriter.write();
        idaPW.close();
        // use the idastandalone solver
        File idaOutputFile = new File(outDir, filePrefix + SimDataConstants.IDA_DATA_EXTENSION);
        // String sundialsSolverExecutable = "C:\\Developer\\Eclipse\\workspace\\VCell 4.8\\SundialsSolverStandalone_NoMessaging.exe";
        String executableName = null;
        try {
            executableName = SolverUtilities.getExes(SolverDescription.IDA)[0].getAbsolutePath();
        } catch (IOException e) {
            throw new RuntimeException("failed to get executable for solver " + SolverDescription.IDA.getDisplayLabel() + ": " + e.getMessage(), e);
        }
        Executable executable = new Executable(new String[] { executableName, idaInputFile.getAbsolutePath(), idaOutputFile.getAbsolutePath() });
        executable.start();
        /*			// Generate .cvodeInput string
			File cvodeFile = new File(outDir,filePrefix+SimDataConstants.CVODEINPUT_DATA_EXTENSION);
			PrintWriter cvodePW = new java.io.PrintWriter(cvodeFile);
			SimulationJob simJob = new SimulationJob(sim, 0, null);
		    CVodeFileWriter cvodeFileWriter = new CVodeFileWriter(cvodePW, simJob);
			cvodeFileWriter.write();
			cvodePW.close();

			// use the cvodeStandalone solver
			File cvodeOutputFile = new File(outDir,filePrefix+SimDataConstants.IDA_DATA_EXTENSION);
			String sundialsSolverExecutable = PropertyLoader.getRequiredProperty(PropertyLoader.sundialsSolverExecutableProperty);
			Executable executable = new Executable(new String[]{sundialsSolverExecutable, cvodeFile.getAbsolutePath(), cvodeOutputFile.getAbsolutePath()});
			executable.start();
*/
        // get the result
        ODESolverResultSet odeSolverResultSet = getODESolverResultSet(simJob, idaOutputFile.getPath());
        // remove CVOde input and output files ??
        idaInputFile.delete();
        idaOutputFile.delete();
        // 
        // print header
        // 
        File outputFile = new File(outDir, "results" + filePrefix + ".csv");
        java.io.PrintStream outputStream = new java.io.PrintStream(new java.io.BufferedOutputStream(new java.io.FileOutputStream(outputFile)));
        outputStream.print("time");
        for (int i = 0; i < testSpec.getVarsList().length; i++) {
            outputStream.print("," + testSpec.getVarsList()[i]);
        }
        outputStream.println();
        // 
        // extract data for time and species
        // 
        double[][] data = new double[testSpec.getVarsList().length + 1][];
        int column = odeSolverResultSet.findColumn("t");
        data[0] = odeSolverResultSet.extractColumn(column);
        int origDataLength = data[0].length;
        for (int i = 0; i < testSpec.getVarsList().length; i++) {
            column = odeSolverResultSet.findColumn(testSpec.getVarsList()[i]);
            if (column == -1) {
                Variable var = simJob.getSimulationSymbolTable().getVariable(testSpec.getVarsList()[i]);
                data[i + 1] = new double[data[0].length];
                if (var instanceof cbit.vcell.math.Constant) {
                    double value = ((cbit.vcell.math.Constant) var).getExpression().evaluateConstant();
                    for (int j = 0; j < data[i + 1].length; j++) {
                        data[i + 1][j] = value;
                    }
                } else {
                    throw new RuntimeException("Did not find " + testSpec.getVarsList()[i] + " in simulation");
                }
            } else {
                data[i + 1] = odeSolverResultSet.extractColumn(column);
            }
        }
        // 
        // for each time, print row
        // 
        int index = 0;
        double[] sampleTimes = new double[testSpec.getNumTimeSteps() + 1];
        for (int i = 0; i <= testSpec.getNumTimeSteps(); i++) {
            sampleTimes[i] = endTime * i / testSpec.getNumTimeSteps();
        }
        Model vcModel = bioModel.getModel();
        ReservedSymbol kMole = vcModel.getKMOLE();
        for (int i = 0; i < sampleTimes.length; i++) {
            // 
            while (true) {
                // 
                if (index == odeSolverResultSet.getRowCount() - 1) {
                    if (data[0][index] == sampleTimes[i]) {
                        break;
                    } else {
                        throw new RuntimeException("sampleTime does not match at last time point");
                    }
                }
                // 
                if (data[0][index + 1] > sampleTimes[i]) {
                    break;
                }
                // 
                // sampleTime must be later in our data list.
                // 
                index++;
            }
            // if data[0][index] == sampleTime no need to interpolate
            if (data[0][index] == sampleTimes[i]) {
                // if timeFactor is not 1.0, time is not in seconds (mins or hrs); if timeFactor is 60, divide sampleTime/60; if it is 3600, divide sampleTime/3600.
                // if (timeFactor != 1.0) {
                // outputStream.print(data[0][index]/timeFactor);
                // } else {
                outputStream.print(data[0][index]);
                // }
                for (int j = 0; j < testSpec.getVarsList().length; j++) {
                    // SBMLImporter.SBVCConcentrationUnits spConcUnits = speciesUnitsHash.get(testSpec.getVarsList()[j]);
                    // if (spConcUnits != null) {
                    // VCUnitDefinition sbunits = spConcUnits.getSBConcentrationUnits();
                    // VCUnitDefinition vcunits = spConcUnits.getVCConcentrationUnits();
                    // SBMLUnitParameter unitFactor = SBMLUtils.getConcUnitFactor("spConcParam", vcunits, sbunits, kMole);
                    // outputStream.print("," + data[j + 1][index] * unitFactor.getExpression().evaluateConstant()); 		//earlier, hack unitfactor = 0.000001
                    // earlier, hack unitfactor = 0.000001
                    outputStream.print("," + data[j + 1][index]);
                // }
                }
                // System.out.println("No interpolation needed!");
                outputStream.println();
            } else {
                // if data[0][index] < sampleTime, must interpolate
                double fraction = (sampleTimes[i] - data[0][index]) / (data[0][index + 1] - data[0][index]);
                // if timeFactor is not 1.0, time is not in seconds (mins or hrs); if timeFactor is 60, divide sampleTime/60; if it is 3600, divide sampleTime/3600.
                // if (timeFactor != 1.0) {
                // outputStream.print(sampleTimes[i]/timeFactor);
                // } else {
                outputStream.print(sampleTimes[i]);
                // }
                for (int j = 0; j < testSpec.getVarsList().length; j++) {
                    double interpolatedValue = 0.0;
                    double[] speciesVals = null;
                    double[] times = null;
                    // Currently using 2nd order interpolation
                    if (index == 0) {
                        // can only do 1st order interpolation
                        times = new double[] { data[0][index], data[0][index + 1] };
                        speciesVals = new double[] { data[j + 1][index], data[j + 1][index + 1] };
                        interpolatedValue = MathTestingUtilities.taylorInterpolation(sampleTimes[i], times, speciesVals);
                    } else if (index >= 1 && index <= origDataLength - 3) {
                        double val_1 = Math.abs(sampleTimes[i] - data[0][index - 1]);
                        double val_2 = Math.abs(sampleTimes[i] - data[0][index + 2]);
                        if (val_1 < val_2) {
                            times = new double[] { data[0][index - 1], data[0][index], data[0][index + 1] };
                            speciesVals = new double[] { data[j + 1][index - 1], data[j + 1][index], data[j + 1][index + 1] };
                        } else {
                            times = new double[] { data[0][index], data[0][index + 1], data[0][index + 2] };
                            speciesVals = new double[] { data[j + 1][index], data[j + 1][index + 1], data[j + 1][index + 2] };
                        }
                        interpolatedValue = MathTestingUtilities.taylorInterpolation(sampleTimes[i], times, speciesVals);
                    } else {
                        times = new double[] { data[0][index - 1], data[0][index], data[0][index + 1] };
                        speciesVals = new double[] { data[j + 1][index - 1], data[j + 1][index], data[j + 1][index + 1] };
                        interpolatedValue = MathTestingUtilities.taylorInterpolation(sampleTimes[i], times, speciesVals);
                    }
                    // // Currently using 1st order interpolation
                    // times = new double[] { data[0][index], data[0][index+1] };
                    // speciesVals = new double[] { data[j+1][index], data[j+1][index+1] };
                    // interpolatedValue = taylorInterpolation(sampleTimes[i], times, speciesVals);
                    // interpolatedValue = interpolatedValue * unitFactor.getExpression().evaluateConstant(); 		//earlier, hack unitfactor = 0.000001
                    // System.out.println("Sample time: " + sampleTimes[i] + ", between time[" + index + "]=" + data[0][index]+" and time["+(index+1)+"]="+(data[0][index+1])+", interpolated = "+interpolatedValue);
                    outputStream.print("," + interpolatedValue);
                }
                outputStream.println();
            }
        }
        outputStream.close();
        return outputFile;
    } catch (Exception e) {
        e.printStackTrace(System.out);
        // File outputFile = new File(outDir,"results" + filePrefix + ".csv");
        throw new SolverException(e.getMessage());
    }
}

Example 7 with SBMLImporter

use of org.vcell.sbml.vcell.SBMLImporter in project vcell by virtualcell.

the class SBMLImporterTest method main.

public static void main(String[] args) {
    Logging.init();
    try {
        VCLogger vcl = new TLogger();
        SBMLImporter imp = new SBMLImporter("samp.sbml", vcl, false);
        for (; ; ) {
            /*
			JFileChooser jfc = new JFileChooser( new File(System.getProperty("user.dir")));
			int returnVal = jfc.showOpenDialog(null);
			if(returnVal == JFileChooser.APPROVE_OPTION) {
				File f= jfc.getSelectedFile();
				BioModel bm = sa.importSBML(f);
				System.out.println(bm.getName());
			}
			*/
            BioModel bm = imp.getBioModel();
            System.out.println(bm.getName());
        }
    } catch (Exception e) {
        System.err.println(e.getMessage());
        e.printStackTrace();
    }
}

Example 8 with SBMLImporter

use of org.vcell.sbml.vcell.SBMLImporter in project vcell by virtualcell.

the class SBMLSpatialTest method test.

// @Test
public void test() throws Exception {
    // BioModel bioModel1 = BioModelTest.getExampleWithImage();
    URL vcmlURL = SBMLSpatialTest.class.getResource("Solver_Suite_6_2.vcml");
    File vcmlFile = new File(vcmlURL.toURI());
    BioModel bioModel1 = XmlHelper.XMLToBioModel(new XMLSource(vcmlFile));
    bioModel1.refreshDependencies();
    // for (int i = 0; i<bioModel1.getNumSimulationContexts(); i++){
    for (int i = 5; i == 5; i++) {
        SimulationContext sc1 = bioModel1.getSimulationContext(i);
        if (sc1.getApplicationType() != Application.NETWORK_DETERMINISTIC) {
            System.err.println(sc1.getName() + " is not a network determistic application");
            continue;
        }
        boolean isSpatial = sc1.getGeometry().getDimension() > 0;
        SBMLExporter exporter = new SBMLExporter(bioModel1, 3, 1, isSpatial);
        sc1.refreshMathDescription(null, NetworkGenerationRequirements.ComputeFullNoTimeout);
        // sc1.setMathDescription(sc1.createNewMathMapping(null, NetworkGenerationRequirements.ComputeFullNoTimeout).getMathDescription());
        exporter.setSelectedSimContext(sc1);
        VCellSBMLDoc sbmlDoc = exporter.convertToSBML();
        for (UnitDefinition unitDefn : sbmlDoc.model.getListOfUnitDefinitions()) {
            for (Unit unit : unitDefn.getListOfUnits()) {
                System.out.println(unit.getKind());
                if (!unit.isSetKind()) {
                    throw new RuntimeException("kind of unit " + unit.printUnit() + " of UnitDefn " + UnitDefinition.printUnits(unitDefn) + " is not set");
                }
            }
        }
        // sbmlDoc.document.setConsistencyChecks(CHECK_CATEGORY.UNITS_CONSISTENCY, false);
        // int numErrors = sbmlDoc.document.checkConsistency();
        // System.out.println("consistency check, num errors = "+numErrors);
        // if (numErrors>0){
        // SBMLErrorLog errorLog = sbmlDoc.document.getListOfErrors();
        // for (int err=0; err<errorLog.getErrorCount(); err++){
        // System.err.println("ERROR IN EXPORTED SBML: "+errorLog.getError(err).getMessage());
        // }
        // //Assert.fail("generated SBML document was found to be inconsistent");
        // }
        String sbmlString = sbmlDoc.xmlString;
        File tempFile = File.createTempFile("sbmlSpatialTest_SBML_", ".sbml.xml");
        FileUtils.write(tempFile, sbmlString);
        System.out.println(tempFile);
        try {
            VCLogger argVCLogger = new TLogger();
            SBMLImporter importer = new SBMLImporter(tempFile.getAbsolutePath(), argVCLogger, isSpatial);
            BioModel bioModel2 = importer.getBioModel();
            File tempFile2 = File.createTempFile("sbmlSpatialTest_Biomodel_", ".vcml.xml");
            FileUtils.write(tempFile2, XmlHelper.bioModelToXML(bioModel2));
            System.out.println(tempFile2);
            // if (true) { throw new RuntimeException("stop"); }
            bioModel2.refreshDependencies();
            SimulationContext sc2 = bioModel2.getSimulationContext(0);
            // sc2.refreshMathDescription(null, NetworkGenerationRequirements.ComputeFullNoTimeout);
            sc2.setMathDescription(sc2.createNewMathMapping(null, NetworkGenerationRequirements.ComputeFullNoTimeout).getMathDescription());
            if (!sc1.getMathDescription().isValid()) {
                throw new RuntimeException("sc1.math is not valid");
            }
            if (!sc2.getMathDescription().isValid()) {
                throw new RuntimeException("sc2.math is not valid");
            }
            MathCompareResults mathCompareResults = MathDescription.testEquivalency(SimulationSymbolTable.createMathSymbolTableFactory(), sc1.getMathDescription(), sc2.getMathDescription());
            if (!mathCompareResults.isEquivalent()) {
                System.out.println("MATH DESCRIPTION 1 <UNCHANGED>");
                System.out.println(sc1.getMathDescription().getVCML_database());
                System.out.println("MATH DESCRIPTION 2 <UNCHANGED>");
                System.out.println(sc2.getMathDescription().getVCML_database());
                // if (mathCompareResults.decision  == Decision.MathDifferent_SUBDOMAINS_DONT_MATCH){
                // BioModel bioModel1_copy = XmlHelper.XMLToBioModel(new XMLSource(vcmlFile));
                // bioModel1_copy.refreshDependencies();
                // SimulationContext sc1_copy = bioModel1_copy.getSimulationContext(i);
                // VCImage image = sc1_copy.getGeometry().getGeometrySpec().getImage();
                // if (image!=null){
                // ArrayList<VCPixelClass> pcList = new ArrayList<VCPixelClass>();
                // for (VCPixelClass pc : image.getPixelClasses()){
                // pcList.add(new VCPixelClass(pc.getKey(),SBMLExporter.DOMAIN_TYPE_PREFIX+pc.getPixelClassName(),pc.getPixel()));
                // }
                // image.setPixelClasses(pcList.toArray(new VCPixelClass[0]));
                // }
                // for (GeometryClass gc : sc1_copy.getGeometry().getGeometryClasses()){
                // System.out.println("name before "+gc.getName());
                // gc.setName(SBMLExporter.DOMAIN_TYPE_PREFIX+gc.getName());
                // System.out.println("name after "+gc.getName());
                // }
                // sc1_copy.checkValidity();
                // bioModel1_copy.refreshDependencies();
                // sc1_copy.getGeometry().precomputeAll(new GeometryThumbnailImageFactoryAWT(), true, true);
                // sc1_copy.setMathDescription(sc1_copy.createNewMathMapping(null, NetworkGenerationRequirements.ComputeFullNoTimeout).getMathDescription());
                // MathCompareResults mathCompareResults_renamedSubdomains = MathDescription.testEquivalency(SimulationSymbolTable.createMathSymbolTableFactory(),sc1_copy.getMathDescription(), sc2.getMathDescription());
                // if (!mathCompareResults_renamedSubdomains.isEquivalent()){
                // System.out.println("MATH DESCRIPTION 1 <RENAMED>");
                // System.out.println(sc1_copy.getMathDescription().getVCML_database());
                // Assert.fail(mathCompareResults_renamedSubdomains.decision+" "+mathCompareResults_renamedSubdomains.details);
                // }
                // }else{
                System.err.println(mathCompareResults.decision + " " + mathCompareResults.details);
            // }
            } else {
                System.out.println("MATHS WERE EQUIVALENT");
            }
        } finally {
            tempFile.delete();
        }
    }
    // loop over determinstic applications
    System.out.println("done");
}