diff --git a/moose-core/Makefile b/moose-core/Makefile index 65b7c9b361a365aaac18ec489b87d74ac9c718ac..5dda459827e9ab9c719b762b9e19946bd0ab0daf 100644 --- a/moose-core/Makefile +++ b/moose-core/Makefile @@ -52,15 +52,11 @@ # generating the python interface using swig. # # USE_MPI - compile with support for parallel computing through MPICH library -# -# USE_SBML (default value: 0) - compile with support for the Systems Biology -# Markup Language (SBML). This allows you to read and write chemical -# kinetic models in the simulator-indpendent SBML format. -# + # Default values for flags. The operator ?= assigns the given value only if the # variable is not already defined. -#USE_SBML?=0 + USE_HDF5?=1 USE_CUDA?=0 USE_NEUROKIT?=0 diff --git a/moose-core/python/moose/SBML/readSBML.py b/moose-core/python/moose/SBML/readSBML.py index 2765c5e045a4d600aea8a23c0398e422c8bd5c96..801d3681ec23536b60dce134f2faa67480c70571 100644 --- a/moose-core/python/moose/SBML/readSBML.py +++ b/moose-core/python/moose/SBML/readSBML.py @@ -10,10 +10,10 @@ ** This program is part of 'MOOSE', the ** Messaging Object Oriented Simulation Environment, ** also known as GENESIS 3 base code. -** copyright (C) 2003-2016 Upinder S. Bhalla. and NCBS +** copyright (C) 2003-2017 Upinder S. Bhalla. and NCBS Created : Thu May 12 10:19:00 2016(+0530) Version -Last-Updated: Wed Sep 28 +Last-Updated: Wed Jan 11 2017 By: **********************************************************************/ @@ -23,7 +23,7 @@ import sys import os.path import collections import moose - +from validation import validateModel ''' TODO in @@ -57,7 +57,6 @@ try: except ImportError: pass - def mooseReadSBML(filepath, loadpath, solver="ee"): global foundLibSBML_ if not foundLibSBML_: @@ -74,12 +73,9 @@ def mooseReadSBML(filepath, loadpath, solver="ee"): with open(filepath, "r") as filep: filep = open(filepath, "r") document = libsbml.readSBML(filepath) - num_errors = document.getNumErrors() - if (num_errors > 0): - print("Encountered the following SBML errors:") - document.printErrors() - return moose.element('/') - else: + tobecontinue = False + tobecontinue = validateModel(document) + if tobecontinue: level = document.getLevel() version = document.getVersion() print(("\n" + "File: " + filepath + " (Level " + @@ -129,8 +125,10 @@ def mooseReadSBML(filepath, loadpath, solver="ee"): global comptSbmlidMooseIdMap global warning warning = " " + global msg + msg = " " comptSbmlidMooseIdMap = {} - print(("modelPath:" + basePath.path)) + #print(("modelPath:" + basePath.path)) globparameterIdValue = {} modelAnnotaInfo = {} mapParameter(model, globparameterIdValue) @@ -138,8 +136,9 @@ def mooseReadSBML(filepath, loadpath, solver="ee"): basePath, model, comptSbmlidMooseIdMap) if errorFlag: specInfoMap = {} - errorFlag = createSpecies( + errorFlag,warning = createSpecies( basePath, model, comptSbmlidMooseIdMap, specInfoMap, modelAnnotaInfo) + #print(errorFlag,warning) if errorFlag: errorFlag = createRules( model, specInfoMap, globparameterIdValue) @@ -148,7 +147,6 @@ def mooseReadSBML(filepath, loadpath, solver="ee"): model, specInfoMap, modelAnnotaInfo, globparameterIdValue) getModelAnnotation( model, baseId, basePath) - if not errorFlag: print(msg) # Any time in the middle if SBML does not read then I @@ -157,7 +155,11 @@ def mooseReadSBML(filepath, loadpath, solver="ee"): # built which is not correct print "Deleted rest of the # model" moose.delete(basePath) - return baseId + basePath = moose.Shell('/') + return basePath + else: + print("Validation failed while reading the model.") + return moose.element('/') def setupEnzymaticReaction(enz, groupName, enzName, @@ -509,9 +511,7 @@ def createReaction(model, specInfoMap, modelAnnotaInfo, globparameterIdValue): nummodifiers = reac.getNumModifiers() if not (numRcts and numPdts): - print( - rName, - " : Substrate and Product is missing, we will be skiping creating this reaction in MOOSE") + print("Warning: %s" %(rName)," : Substrate or Product is missing, we will be skiping creating this reaction in MOOSE") reactionCreated = False elif (reac.getNumModifiers() > 0): reactionCreated, reaction_ = setupMMEnzymeReaction( @@ -805,9 +805,8 @@ def createSpecies(basePath, model, comptSbmlidMooseIdMap, # - Need to add group name if exist in pool # - Notes # print "species " - if not (model.getNumSpecies()): - return False + return (False,"number of species is zero") else: for sindex in range(0, model.getNumSpecies()): spe = model.getSpecies(sindex) @@ -876,8 +875,7 @@ def createSpecies(basePath, model, comptSbmlidMooseIdMap, initvalue = initvalue * unitfactor elif spe.isSetInitialConcentration(): initvalue = spe.getInitialConcentration() - print( - " Since hasonlySubUnit is true and concentration is set units are not checked") + print(" Since hasonlySubUnit is true and concentration is set units are not checked") poolId.nInit = initvalue elif hasonlySubUnit == False: @@ -885,8 +883,7 @@ def createSpecies(basePath, model, comptSbmlidMooseIdMap, if spe.isSetInitialAmount(): initvalue = spe.getInitialAmount() # initAmount is set we need to convert to concentration - initvalue = initvalue / \ - comptSbmlidMooseIdMap[comptId]["size"] + initvalue = initvalue / comptSbmlidMooseIdMap[comptId]["size"] elif spe.isSetInitialConcentration(): initvalue = spe.getInitialConcentration() @@ -911,9 +908,9 @@ def createSpecies(basePath, model, comptSbmlidMooseIdMap, print( "Invalid SBML: Either initialConcentration or initialAmount must be set or it should be found in assignmentRule but non happening for ", sName) - return False + return (False,"Invalid SBML: Either initialConcentration or initialAmount must be set or it should be found in assignmentRule but non happening for ",sName) - return True + return (True," ") def transformUnit(unitForObject, hasonlySubUnit=False): @@ -995,7 +992,7 @@ def createCompartment(basePath, model, comptSbmlidMooseIdMap): # ToDoList : Check what should be done for the spaitialdimension is 2 or # 1, area or length if not(model.getNumCompartments()): - return False + return False, else: for c in range(0, model.getNumCompartments()): compt = model.getCompartment(c) @@ -1115,19 +1112,25 @@ def findCompartment(element): return element if __name__ == "__main__": - - filepath = sys.argv[1] - path = sys.argv[2] - - f = open(filepath, 'r') - - if path == '': - loadpath = filepath[filepath.rfind('/'):filepath.find('.')] + try: + sys.argv[1] + except IndexError: + print("Filename or path not given") + exit(0) else: - loadpath = path - - read = mooseReadSBML(filepath, loadpath) - if read: - print(" Read to path", loadpath) - else: - print(" could not read SBML to MOOSE") + filepath = sys.argv[1] + if not os.path.exists(filepath): + print("Filename or path does not exist",filepath) + + else: + try: + sys.argv[2] + except : + modelpath = filepath[filepath.rfind('/'):filepath.find('.')] + else: + modelpath = sys.argv[2] + read = mooseReadSBML(filepath, modelpath) + if read: + print(" Model read to moose path "+ modelpath) + else: + print(" could not read SBML to MOOSE") \ No newline at end of file diff --git a/moose-core/python/moose/SBML/validation.py b/moose-core/python/moose/SBML/validation.py new file mode 100644 index 0000000000000000000000000000000000000000..41ab4140b6a1a1ae659d2738fd0471e793756fcb --- /dev/null +++ b/moose-core/python/moose/SBML/validation.py @@ -0,0 +1,95 @@ +foundLibSBML_ = False +try: + from libsbml import * + foundLibSBML_ = True +except Exception as e: + pass + +def validateModel(sbmlDoc): + if sbmlDoc.getNumErrors() > 0: + tobecontinued = False + for i in range(0,sbmlDoc.getNumErrors()): + print (sbmlDoc.getError(i).getMessage()) + return False + + if (not sbmlDoc): + print("validateModel: given a null SBML Document") + return False + + consistencyMessages = "" + validationMessages = "" + noProblems = True + numCheckFailures = 0 + numConsistencyErrors = 0 + numConsistencyWarnings = 0 + numValidationErrors = 0 + numValidationWarnings = 0 + # Once the whole model is done and before it gets written out, + # it's important to check that the whole model is in fact complete, + # consistent and valid. + numCheckFailures = sbmlDoc.checkInternalConsistency() + if (numCheckFailures > 0): + for i in range(0, numCheckFailures): + sbmlErr = sbmlDoc.getError(i) + if (sbmlErr.isFatal() or sbmlErr.isError()): + noProblems = False + numConsistencyErrors += 1 + else: + numConsistencyWarnings += 1 + constStr = sbmlDoc.printErrors() + if sbmlDoc.printErrors(): + consistencyMessages = constStr + + # If the internal checks fail, it makes little sense to attempt + # further validation, because the model may be too compromised to + # be properly interpreted. + + if (numConsistencyErrors > 0): + consistencyMessages += "Further validation aborted." + else: + numCheckFailures = sbmlDoc.checkConsistency() + validationMessages; + #numCheckFailures = sbmlDoc.checkL3v1Compatibility() + if (numCheckFailures > 0): + for i in range(0, (numCheckFailures)): + consistencyMessages = sbmlDoc.getErrorLog().toString() + sbmlErr = sbmlDoc.getError(i) + if (sbmlErr.isFatal() or sbmlErr.isError()): + noProblems = False + numValidationErrors += 1 + else: + numValidationWarnings += 1 + + warning = sbmlDoc.getErrorLog().toString() + oss = sbmlDoc.printErrors() + validationMessages = oss + + if (noProblems): + return True + else: + if consistencyMessages is None: + consistencyMessages = "" + if consistencyMessages != "": + print("consistency Warning: " + consistencyMessages) + + if (numConsistencyErrors > 0): + print("ERROR: encountered " +str(numConsistencyErrors) +" consistency error in model " +sbmlDoc.getModel().getId() +"'.") + + if (numConsistencyWarnings > 0): + print("Notice: encountered " +str(numConsistencyWarnings) +" consistency warning in model " +sbmlDoc.getModel().getId() +"'.") + + if (numValidationErrors > 0): + print("ERROR: encountered " + str(numValidationErrors) + " validation error in model " +sbmlDoc.getModel().getId() + "'.") + if (numValidationWarnings > 0): + print("Notice: encountered " +str(numValidationWarnings) +" validation warning in model " +sbmlDoc.getModel().getId() +"'.") + + if validationMessages: + print(validationMessages) + + return False + # return ( numConsistencyErrors == 0 and numValidationErrors == 0, + # consistencyMessages) + +if __name__ == '__main__': + sbmlDoc = libsbml.readSBML('00001-sbml-l3v1.xml') + validateModel(sbmlDoc) diff --git a/moose-core/python/moose/SBML/writeSBML.py b/moose-core/python/moose/SBML/writeSBML.py index 3057f1002b3d356114f579a195a97c0804f568b1..bcfaf2eb2824b52337fd8b058f7b2dd5ba4d5fc8 100644 --- a/moose-core/python/moose/SBML/writeSBML.py +++ b/moose-core/python/moose/SBML/writeSBML.py @@ -9,26 +9,28 @@ ** This program is part of 'MOOSE', the ** Messaging Object Oriented Simulation Environment, ** also known as GENESIS 3 base code. -** copyright (C) 2003-2016 Upinder S. Bhalla. and NCBS +** copyright (C) 2003-2017 Upinder S. Bhalla. and NCBS Created : Friday May 27 12:19:00 2016(+0530) Version -Last-Updated: Thursday Oct 27 11:20:00 2016(+0530) - By: +Last-Updated: Wed Jan 11 15:20:00 2017(+0530) + By: **********************************************************************/ /**************************** ''' - -import moose +import sys import re from collections import Counter -import networkx as nx -import matplotlib.pyplot as plt -import sys + +import moose +from validation import validateModel +from moose.chemUtil.chemConnectUtil import * +from moose.chemUtil.graphUtils import * + # ToDo: -# Table should be written -# Group's should be added +# Table should be written +# Group's should be added # boundary condition for buffer pool having assignment statment constant # shd be false @@ -41,10 +43,11 @@ except Exception as e: def mooseWriteSBML(modelpath, filename, sceneitems={}): global foundLibSBML_ + msg = " " if not foundLibSBML_: print('No python-libsbml found.' '\nThis module can be installed by following command in terminal:' - '\n\t easy_install python-libsbl' + '\n\t easy_install python-libsbml' ) return -1, msg @@ -57,12 +60,11 @@ def mooseWriteSBML(modelpath, filename, sceneitems={}): # validatemodel sbmlOk = False - global spe_constTrue, cmin, cmax + global spe_constTrue spe_constTrue = [] global nameList_ nameList_ = [] - - autoCoordinateslayout = False + positionInfoexist = False xmlns = XMLNamespaces() xmlns.add("http://www.sbml.org/sbml/level3/version1") xmlns.add("http://www.moose.ncbs.res.in", "moose") @@ -74,55 +76,57 @@ def mooseWriteSBML(modelpath, filename, sceneitems={}): cremodel_.setExtentUnits("substance") cremodel_.setSubstanceUnits("substance") neutralNotes = "" + specieslist = moose.wildcardFind(modelpath + '/##[ISA=PoolBase]') - neutralPath = getGroupinfo(specieslist[0]) - if moose.exists(neutralPath.path + '/info'): - neutralInfo = moose.element(neutralPath.path + '/info') - neutralNotes = neutralInfo.notes - if neutralNotes != "": - cleanNotes = convertNotesSpecialChar(neutralNotes) - notesString = "<body xmlns=\"http://www.w3.org/1999/xhtml\">\n \t \t" + \ - neutralNotes + "\n\t </body>" - cremodel_.setNotes(notesString) - srcdesConnection = {} - cmin, cmax = 0, 1 + if specieslist: + neutralPath = getGroupinfo(specieslist[0]) + if moose.exists(neutralPath.path + '/info'): + neutralInfo = moose.element(neutralPath.path + '/info') + neutralNotes = neutralInfo.notes + if neutralNotes != "": + cleanNotes = convertNotesSpecialChar(neutralNotes) + notesString = "<body xmlns=\"http://www.w3.org/1999/xhtml\">\n \t \t" + \ + neutralNotes + "\n\t </body>" + cremodel_.setNotes(notesString) + srcdesConnection = {} if not bool(sceneitems): - autoCoordinateslayout = True - srcdesConnection = setupItem(modelpath) - meshEntry = setupMeshObj(modelpath) - cmin, cmax, sceneitems = autoCoordinates( - meshEntry, srcdesConnection) - + meshEntry,xmin,xmax,ymin,ymax,positionInfoexist,sitem = setupMeshObj(modelpath) + #moose.coordinateUtil.setupItem(modelpath,srcdesConnection) + setupItem(modelpath,srcdesConnection) + if not positionInfoexist: + autoCoordinates(meshEntry,srcdesConnection) + writeUnits(cremodel_) modelAnno = writeSimulationAnnotation(modelpath) if modelAnno: cremodel_.setAnnotation(modelAnno) + compartexist = writeCompt(modelpath, cremodel_) - species = writeSpecies( - modelpath, - cremodel_, - sbmlDoc, - sceneitems, - autoCoordinateslayout) - if species: - writeFunc(modelpath, cremodel_) - writeReac(modelpath, cremodel_, sceneitems, autoCoordinateslayout) - writeEnz(modelpath, cremodel_, sceneitems, autoCoordinateslayout) - - consistencyMessages = "" - SBMLok = validateModel(sbmlDoc) - if (SBMLok): - writeTofile = filepath + "/" + filename + '.xml' - writeSBMLToFile(sbmlDoc, writeTofile) - return True, consistencyMessages, writeTofile - - if (not SBMLok): - cerr << "Errors encountered " << endl - return -1, consistencyMessages - - -def writeEnz(modelpath, cremodel_, sceneitems, autoCoordinateslayout): + if compartexist == True: + species = writeSpecies( modelpath,cremodel_,sbmlDoc,sceneitems) + if species: + writeFunc(modelpath, cremodel_) + writeReac(modelpath, cremodel_, sceneitems) + writeEnz(modelpath, cremodel_, sceneitems) + consistencyMessages = "" + SBMLok = validateModel(sbmlDoc) + if (SBMLok): + writeTofile = filepath + "/" + filename + '.xml' + writeSBMLToFile(sbmlDoc, writeTofile) + return True, consistencyMessages, writeTofile + + if (not SBMLok): + cerr << "Errors encountered " << endl + return -1, consistencyMessages + else: + return False, "Atleast one compartment should exist to write SBML" + +def calPrime(x): + prime = int((20 * (float(x - cmin) / float(cmax - cmin))) - 10) + return prime + +def writeEnz(modelpath, cremodel_, sceneitems): for enz in moose.wildcardFind(modelpath + '/##[ISA=EnzBase]'): enzannoexist = False enzGpnCorCol = " " @@ -135,25 +139,26 @@ def writeEnz(modelpath, cremodel_, sceneitems, autoCoordinateslayout): notesE = Anno.notes element = moose.element(enz) ele = getGroupinfo(element) - if ele.className == "Neutral" or sceneitems[ - element] or Anno.color or Anno.textColor: - enzannoexist = True - + if element.className == "Neutral" or sceneitems or Anno.x or Anno.y: + enzannoexist = True if enzannoexist: + enzAnno = "<moose:ModelAnnotation>\n" if ele.className == "Neutral": - enzGpnCorCol = "<moose:Group> " + ele.name + " </moose:Group>\n" - if sceneitems[element]: - v = sceneitems[element] - if autoCoordinateslayout == False: - enzGpnCorCol = enzGpnCorCol + "<moose:xCord>" + \ - str(v['x']) + "</moose:xCord>\n" + \ - "<moose:yCord>" + str(v['y']) + "</moose:yCord>\n" - elif autoCoordinateslayout == True: - x = calPrime(v['x']) - y = calPrime(v['y']) - enzGpnCorCol = enzGpnCorCol + "<moose:xCord>" + \ - str(x) + "</moose:xCord>\n" + \ - "<moose:yCord>" + str(y) + "</moose:yCord>\n" + enzGpnCorCol = "<moose:Group>" + ele.name + "</moose:Group>\n" + if sceneitems: + #Saved from GUI, then scene co-ordinates are passed + enzGpnCorCol = enzGpnCorCol + "<moose:xCord>" + \ + str(sceneitems[enz]['x']) + "</moose:xCord>\n" + \ + "<moose:yCord>" + \ + str(sceneitems[enz]['y'])+ "</moose:yCord>\n" + else: + #Saved from cmdline,genesis coordinates are kept as its + # SBML, cspace, python, then auto-coordinates are done + #and coordinates are updated in moose Annotation field + enzGpnCorCol = enzGpnCorCol + "<moose:xCord>" + \ + str(Anno.x) + "</moose:xCord>\n" + \ + "<moose:yCord>" + \ + str(Anno.y)+ "</moose:yCord>\n" if Anno.color: enzGpnCorCol = enzGpnCorCol + "<moose:bgColor>" + Anno.color + "</moose:bgColor>\n" if Anno.textColor: @@ -501,7 +506,7 @@ def listofname(reacSub, mobjEnz): nameList_.append(clean_name) -def writeReac(modelpath, cremodel_, sceneitems, autoCoordinateslayout): +def writeReac(modelpath, cremodel_, sceneitems): for reac in moose.wildcardFind(modelpath + '/##[ISA=ReacBase]'): reacSub = reac.neighbors["sub"] reacPrd = reac.neighbors["prd"] @@ -537,34 +542,34 @@ def writeReac(modelpath, cremodel_, sceneitems, autoCoordinateslayout): reaction.setNotes(notesStringR) element = moose.element(reac) ele = getGroupinfo(element) - if ele.className == "Neutral" or sceneitems[ - element] or Anno.color or Anno.textColor: + if element.className == "Neutral" or sceneitems or Anno.x or Anno.y: reacannoexist = True if reacannoexist: reacAnno = "<moose:ModelAnnotation>\n" if ele.className == "Neutral": reacAnno = reacAnno + "<moose:Group>" + ele.name + "</moose:Group>\n" - if sceneitems[element]: - v = sceneitems[element] - if autoCoordinateslayout == False: - reacAnno = reacAnno + "<moose:xCord>" + \ - str(v['x']) + "</moose:xCord>\n" + \ + if sceneitems: + #Saved from GUI, then scene co-ordinates are passed + reacAnno = reacAnno + "<moose:xCord>" + \ + str(sceneitems[reac]['x']) + "</moose:xCord>\n" + \ + "<moose:yCord>" + \ + str(sceneitems[reac]['y'])+ "</moose:yCord>\n" + else: + #Saved from cmdline,genesis coordinates are kept as its + # SBML, cspace, python, then auto-coordinates are done + #and coordinates are updated in moose Annotation field + reacAnno = reacAnno + "<moose:xCord>" + \ + str(Anno.x) + "</moose:xCord>\n" + \ "<moose:yCord>" + \ - str(v['y']) + "</moose:yCord>\n" - elif autoCoordinateslayout == True: - x = calPrime(v['x']) - y = calPrime(v['y']) - reacAnno = reacAnno + "<moose:xCord>" + \ - str(x) + "</moose:xCord>\n" + \ - "<moose:yCord>" + str(y) + "</moose:yCord>\n" + str(Anno.y)+ "</moose:yCord>\n" if Anno.color: reacAnno = reacAnno + "<moose:bgColor>" + Anno.color + "</moose:bgColor>\n" if Anno.textColor: reacAnno = reacAnno + "<moose:textColor>" + \ Anno.textColor + "</moose:textColor>\n" reacAnno = reacAnno + "</moose:ModelAnnotation>" - # s1.appendAnnotation(XMLNode.convertStringToXMLNode(speciAnno)) reaction.setAnnotation(reacAnno) + kl_s, sRL, pRL, compt = "", "", "", "" if not reacSub and not reacPrd: @@ -712,8 +717,7 @@ def convertSpecialChar(str1): return str1 -def writeSpecies(modelpath, cremodel_, sbmlDoc, - sceneitems, autoCoordinateslayout): +def writeSpecies(modelpath, cremodel_, sbmlDoc, sceneitems): # getting all the species for spe in moose.wildcardFind(modelpath + '/##[ISA=PoolBase]'): sName = convertSpecialChar(spe.name) @@ -785,28 +789,29 @@ def writeSpecies(modelpath, cremodel_, sbmlDoc, cleanNotesS + "\n\t </body>" s1.setNotes(notesStringS) + element = moose.element(spe) ele = getGroupinfo(element) - if ele.className == "Neutral" or sceneitems[ - element] or Anno.color or Anno.textColor: + if element.className == "Neutral" or Anno.color or Anno.textColor or sceneitems or Anno.x or Anno.y: speciannoexist = True if speciannoexist: speciAnno = "<moose:ModelAnnotation>\n" if ele.className == "Neutral": speciAnno = speciAnno + "<moose:Group>" + ele.name + "</moose:Group>\n" - if sceneitems[element]: - v = sceneitems[element] - if autoCoordinateslayout == False: - speciAnno = speciAnno + "<moose:xCord>" + \ - str(v['x']) + "</moose:xCord>\n" + \ + if sceneitems: + #Saved from GUI, then scene co-ordinates are passed + speciAnno = speciAnno + "<moose:xCord>" + \ + str(sceneitems[spe]['x']) + "</moose:xCord>\n" + \ + "<moose:yCord>" + \ + str(sceneitems[spe]['y'])+ "</moose:yCord>\n" + else: + #Saved from cmdline,genesis coordinates are kept as its + # SBML, cspace, python, then auto-coordinates are done + #and coordinates are updated in moose Annotation field + speciAnno = speciAnno + "<moose:xCord>" + \ + str(Anno.x) + "</moose:xCord>\n" + \ "<moose:yCord>" + \ - str(v['y']) + "</moose:yCord>\n" - elif autoCoordinateslayout == True: - x = calPrime(v['x']) - y = calPrime(v['y']) - speciAnno = speciAnno + "<moose:xCord>" + \ - str(x) + "</moose:xCord>\n" + \ - "<moose:yCord>" + str(y) + "</moose:yCord>\n" + str(Anno.y)+ "</moose:yCord>\n" if Anno.color: speciAnno = speciAnno + "<moose:bgColor>" + Anno.color + "</moose:bgColor>\n" if Anno.textColor: @@ -819,7 +824,8 @@ def writeSpecies(modelpath, cremodel_, sbmlDoc, def writeCompt(modelpath, cremodel_): # getting all the compartments - for compt in moose.wildcardFind(modelpath + '/##[ISA=ChemCompt]'): + compts = moose.wildcardFind(modelpath + '/##[ISA=ChemCompt]') + for compt in compts: comptName = convertSpecialChar(compt.name) # converting m3 to litre size = compt.volume * pow(10, 3) @@ -836,7 +842,10 @@ def writeCompt(modelpath, cremodel_): c1.setSize(size) c1.setSpatialDimensions(ndim) c1.setUnits('volume') - + if compts: + return True + else: + return False # write Simulation runtime,simdt,plotdt def writeSimulationAnnotation(modelpath): modelAnno = "" @@ -878,7 +887,6 @@ def writeSimulationAnnotation(modelpath): modelAnno = modelAnno + "</moose:ModelAnnotation>" return modelAnno - def writeUnits(cremodel_): unitVol = cremodel_.createUnitDefinition() unitVol.setId("volume") @@ -896,357 +904,28 @@ def writeUnits(cremodel_): unit.setExponent(1.0) unit.setScale(-3) - -def validateModel(sbmlDoc): - # print " sbmlDoc ",sbmlDoc.toSBML() - if (not sbmlDoc): - print("validateModel: given a null SBML Document") - return False - consistencyMessages = "" - validationMessages = "" - noProblems = True - numCheckFailures = 0 - numConsistencyErrors = 0 - numConsistencyWarnings = 0 - numValidationErrors = 0 - numValidationWarnings = 0 - # Once the whole model is done and before it gets written out, - # it's important to check that the whole model is in fact complete, - # consistent and valid. - numCheckFailures = sbmlDoc.checkInternalConsistency() - if (numCheckFailures > 0): - noProblems = False - for i in range(0, numCheckFailures): - sbmlErr = sbmlDoc.getError(i) - if (sbmlErr.isFatal() or sbmlErr.isError()): - numConsistencyErrors += 1 - else: - numConsistencyWarnings += 1 - constStr = sbmlDoc.printErrors() - consistencyMessages = constStr - - # If the internal checks fail, it makes little sense to attempt - # further validation, because the model may be too compromised to - # be properly interpreted. - if (numConsistencyErrors > 0): - consistencyMessages += "Further validation aborted." - else: - numCheckFailures = sbmlDoc.checkConsistency() - #numCheckFailures = sbmlDoc.checkL3v1Compatibility() - if (numCheckFailures > 0): - noProblems = False - for i in range(0, (numCheckFailures)): - consistencyMessages = sbmlDoc.getErrorLog().toString() - sbmlErr = sbmlDoc.getError(i) - if (sbmlErr.isFatal() or sbmlErr.isError()): - numValidationErrors += 1 - else: - numValidationWarnings += 1 - warning = sbmlDoc.getErrorLog().toString() - oss = sbmlDoc.printErrors() - validationMessages = oss - if (noProblems): - return True - else: - if consistencyMessages is None: - consistencyMessages = "" - if consistencyMessages != "": - print(" consistency Warning: " + consistencyMessages) - - if (numConsistencyErrors > 0): - if numConsistencyErrors == 1: - t = "" - else: - t = "s" - print( - "ERROR: encountered " + - numConsistencyErrors + - " consistency error" + - t + - " in model '" + - sbmlDoc.getModel().getId() + - "'.") - if (numConsistencyWarnings > 0): - if numConsistencyWarnings == 1: - t1 = "" - else: - t1 = "s" - print( - "Notice: encountered " + - numConsistencyWarnings + - " consistency warning" + - t + - " in model '" + - sbmlDoc.getModel().getId() + - "'.") - - if (numValidationErrors > 0): - if numValidationErrors == 1: - t2 = "" - else: - t2 = "s" - print( - "ERROR: encountered " + - numValidationErrors + - " validation error" + - t2 + - " in model '" + - sbmlDoc.getModel().getId() + - "'.") - if (numValidationWarnings > 0): - if numValidationWarnings == 1: - t3 = "" - else: - t3 = "s" - - print( - "Notice: encountered " + - numValidationWarnings + - " validation warning" + - t3 + - " in model '" + - sbmlDoc.getModel().getId() + - "'.") - - print(validationMessages) - return (numConsistencyErrors == 0 and numValidationErrors == 0) - # return ( numConsistencyErrors == 0 and numValidationErrors == 0, - # consistencyMessages) - - -def setupItem(modelPath): - '''This function collects information of what is connected to what. \ - eg. substrate and product connectivity to reaction's and enzyme's \ - sumtotal connectivity to its pool are collected ''' - # print " setupItem" - sublist = [] - prdlist = [] - cntDict = {} - zombieType = ['ReacBase', 'EnzBase', 'Function', 'StimulusTable'] - for baseObj in zombieType: - path = '/##[ISA=' + baseObj + ']' - if modelPath != '/': - path = modelPath + path - if ((baseObj == 'ReacBase') or (baseObj == 'EnzBase')): - for items in moose.wildcardFind(path): - sublist = [] - prdlist = [] - uniqItem, countuniqItem = countitems(items, 'subOut') - subNo = uniqItem - for sub in uniqItem: - sublist.append((moose.element(sub), 's', countuniqItem[sub])) - - uniqItem, countuniqItem = countitems(items, 'prd') - prdNo = uniqItem - if (len(subNo) == 0 or len(prdNo) == 0): - print("Substrate Product is empty ", " ", items) - for prd in uniqItem: - prdlist.append((moose.element(prd), 'p', countuniqItem[prd])) - - if (baseObj == 'CplxEnzBase'): - uniqItem, countuniqItem = countitems(items, 'toEnz') - for enzpar in uniqItem: - sublist.append( - (moose.element(enzpar), 't', countuniqItem[enzpar])) - - uniqItem, countuniqItem = countitems(items, 'cplxDest') - for cplx in uniqItem: - prdlist.append( - (moose.element(cplx), 'cplx', countuniqItem[cplx])) - - if (baseObj == 'EnzBase'): - uniqItem, countuniqItem = countitems(items, 'enzDest') - for enzpar in uniqItem: - sublist.append( - (moose.element(enzpar), 't', countuniqItem[enzpar])) - cntDict[items] = sublist, prdlist - elif baseObj == 'Function': - for items in moose.wildcardFind(path): - sublist = [] - prdlist = [] - item = items.path + '/x[0]' - uniqItem, countuniqItem = countitems(item, 'input') - for funcpar in uniqItem: - sublist.append( - (moose.element(funcpar), 'sts', countuniqItem[funcpar])) - - uniqItem, countuniqItem = countitems(items, 'valueOut') - for funcpar in uniqItem: - prdlist.append( - (moose.element(funcpar), 'stp', countuniqItem[funcpar])) - cntDict[items] = sublist, prdlist - else: - for tab in moose.wildcardFind(path): - tablist = [] - uniqItem, countuniqItem = countitems(tab, 'output') - for tabconnect in uniqItem: - tablist.append( - (moose.element(tabconnect), 'tab', countuniqItem[tabconnect])) - cntDict[tab] = tablist - return cntDict - - -def countitems(mitems, objtype): - items = [] - # print "mitems in countitems ",mitems,objtype - items = moose.element(mitems).neighbors[objtype] - uniqItems = set(items) - countuniqItems = Counter(items) - return(uniqItems, countuniqItems) - - -def setupMeshObj(modelRoot): - ''' Setup compartment and its members pool,reaction,enz cplx under self.meshEntry dictionaries \ - self.meshEntry with "key" as compartment, - value is key2:list where key2 represents moose object type,list of objects of a perticular type - e.g self.meshEntry[meshEnt] = { 'reaction': reaction_list,'enzyme':enzyme_list,'pool':poollist,'cplx': cplxlist } - ''' - meshEntry = {} - if meshEntry: - meshEntry.clear() +if __name__ == "__main__": + try: + sys.argv[1] + except IndexError: + print("Filename or path not given") + exit(0) else: - meshEntry = {} - meshEntryWildcard = '/##[ISA=ChemCompt]' - if modelRoot != '/': - meshEntryWildcard = modelRoot + meshEntryWildcard - for meshEnt in moose.wildcardFind(meshEntryWildcard): - mollist = [] - cplxlist = [] - mol_cpl = moose.wildcardFind(meshEnt.path + '/##[ISA=PoolBase]') - funclist = moose.wildcardFind(meshEnt.path + '/##[ISA=Function]') - enzlist = moose.wildcardFind(meshEnt.path + '/##[ISA=EnzBase]') - realist = moose.wildcardFind(meshEnt.path + '/##[ISA=ReacBase]') - tablist = moose.wildcardFind(meshEnt.path + '/##[ISA=StimulusTable]') - if mol_cpl or funclist or enzlist or realist or tablist: - for m in mol_cpl: - if isinstance(moose.element(m.parent), moose.CplxEnzBase): - cplxlist.append(m) - elif isinstance(moose.element(m), moose.PoolBase): - mollist.append(m) - - meshEntry[meshEnt] = {'enzyme': enzlist, - 'reaction': realist, - 'pool': mollist, - 'cplx': cplxlist, - 'table': tablist, - 'function': funclist - } - return(meshEntry) - - -def autoCoordinates(meshEntry, srcdesConnection): - G = nx.Graph() - for cmpt, memb in list(meshEntry.items()): - for enzObj in find_index(memb, 'enzyme'): - G.add_node( - enzObj.path, - label='', - shape='ellipse', - color='', - style='filled', - fontname='Helvetica', - fontsize=12, - fontcolor='blue') - for cmpt, memb in list(meshEntry.items()): - for poolObj in find_index(memb, 'pool'): - #poolinfo = moose.element(poolObj.path+'/info') - G.add_node( - poolObj.path, - label=poolObj.name, - shape='box', - color='', - style='filled', - fontname='Helvetica', - fontsize=12, - fontcolor='blue') - for cplxObj in find_index(memb, 'cplx'): - pass - for reaObj in find_index(memb, 'reaction'): - G.add_node(reaObj.path, label='', shape='record', color='') - - for inn, out in list(srcdesConnection.items()): - if (inn.className == 'ZombieReac'): - arrowcolor = 'green' - elif(inn.className == 'ZombieEnz'): - arrowcolor = 'red' + filepath = sys.argv[1] + if not os.path.exists(filepath): + print("Filename or path does not exist",filepath) else: - arrowcolor = 'blue' - if isinstance(out, tuple): - if len(out[0]) == 0: - print( - inn.className + - ':' + - inn.name + - " doesn't have input message") + try: + sys.argv[2] + except : + modelpath = filepath[filepath.rfind('/'):filepath.find('.')] else: - for items in (items for items in out[0]): - G.add_edge(moose.element(items[0]).path, inn.path) - - if len(out[1]) == 0: - print( - inn.className + - ':' + - inn.name + - "doesn't have output mssg") - else: - for items in (items for items in out[1]): - G.add_edge(inn.path, moose.element(items[0]).path) + modelpath = sys.argv[2] + + moose.loadModel(filepath, modelpath, "gsl") - elif isinstance(out, list): - if len(out) == 0: - print("Func pool doesn't have sumtotal") + written, c, writtentofile = mooseWriteSBML(modelpath, filepath) + if written: + print(" File written to ", writtentofile) else: - for items in (items for items in out): - G.add_edge(moose.element(items[0]).path, inn.path) - #from networkx.drawing.nx_agraph import graphviz_layout - #position = graphviz_layout(G,prog='dot') - - position = nx.pygraphviz_layout(G, prog='dot') - position = nx.spring_layout(G) - #agraph = nx.to_agraph(G) - #agraph.draw("~/out.png", format = 'png', prog = 'dot') - - sceneitems = {} - xycord = [] - cmin = 0 - cmax = 0 - for key, value in list(position.items()): - xycord.append(value[0]) - xycord.append(value[1]) - sceneitems[moose.element(key)] = {'x': value[0], 'y': value[1]} - if len(xycord) > 0: - cmin = min(xycord) - cmax = max(xycord) - return cmin, cmax, sceneitems - - -def calPrime(x): - prime = int((20 * (float(x - cmin) / float(cmax - cmin))) - 10) - return prime - - -def find_index(value, key): - """ Value.get(key) to avoid expection which would raise if empty value in dictionary for a given key """ - if value.get(key) is not None: - return value.get(key) - else: - raise ValueError('no dict with the key found') -if __name__ == "__main__": - - filepath = sys.argv[1] - path = sys.argv[2] - - f = open(filepath, 'r') - - if path == '': - loadpath = filepath[filepath.rfind('/'):filepath.find('.')] - else: - loadpath = path - - moose.loadModel(filepath, loadpath, "gsl") - - written, c, writtentofile = mooseWriteSBML(loadpath, filepath) - if written: - print(" File written to ", writtentofile) - else: - print(" could not write model to SBML file") + print(" could not write model to SBML file") diff --git a/moose-core/python/moose/__init__.py b/moose-core/python/moose/__init__.py index 27c5a3bf26eb4c81d734f870d5c3d0d63e82c889..a1a77c0ea7e188ed7e3ef16b5cf2e246e7861c1c 100644 --- a/moose-core/python/moose/__init__.py +++ b/moose-core/python/moose/__init__.py @@ -27,3 +27,7 @@ from .moose import * # Wrapper to get moose version information. __version__ = moose._moose.__version__ VERSION = moose._moose.VERSION + +import chemUtil.add_Delete_ChemicalSolver +import chemUtil.chemConnectUtil +import chemUtil.graphUtils \ No newline at end of file diff --git a/moose-core/python/moose/chemUtil/__init__.py b/moose-core/python/moose/chemUtil/__init__.py new file mode 100644 index 0000000000000000000000000000000000000000..dce7080f0a0f57a421c2428821689fe824ff50a5 --- /dev/null +++ b/moose-core/python/moose/chemUtil/__init__.py @@ -0,0 +1,3 @@ +from .add_Delete_ChemicalSolver import * +from .chemConnectUtil import * +from .graphUtils import autoCoordinates diff --git a/moose-core/python/moose/add_Delete_ChemicalSolver.py b/moose-core/python/moose/chemUtil/add_Delete_ChemicalSolver.py similarity index 99% rename from moose-core/python/moose/add_Delete_ChemicalSolver.py rename to moose-core/python/moose/chemUtil/add_Delete_ChemicalSolver.py index f788118b7d20ad116c0a3323330ed20a598a712f..6d5886f19bd5a003a5a1fa4a960ad5a5b4774809 100644 --- a/moose-core/python/moose/add_Delete_ChemicalSolver.py +++ b/moose-core/python/moose/chemUtil/add_Delete_ChemicalSolver.py @@ -1,5 +1,4 @@ -from . import _moose as moose - +import moose def moosedeleteChemSolver(modelRoot): """Delete solvers from Chemical Compartment diff --git a/moose-core/python/moose/chemUtil/chemConnectUtil.py b/moose-core/python/moose/chemUtil/chemConnectUtil.py new file mode 100644 index 0000000000000000000000000000000000000000..93d8be9d1442970c5ac748184f6743ba6810fc6b --- /dev/null +++ b/moose-core/python/moose/chemUtil/chemConnectUtil.py @@ -0,0 +1,181 @@ +import moose +import numpy as np +from collections import Counter + +def xyPosition(objInfo,xory): + try: + return(float(moose.element(objInfo).getField(xory))) + except ValueError: + return (float(0)) + +def setupMeshObj(modelRoot): + ''' Setup compartment and its members pool,reaction,enz cplx under self.meshEntry dictionaries \ + self.meshEntry with "key" as compartment, + value is key2:list where key2 represents moose object type,list of objects of a perticular type + e.g self.meshEntry[meshEnt] = { 'reaction': reaction_list,'enzyme':enzyme_list,'pool':poollist,'cplx': cplxlist } + ''' + xmin = 0.0 + xmax = 1.0 + ymin = 0.0 + ymax = 1.0 + listOfitems = {} + positionInfoExist = True + meshEntry = {} + if meshEntry: + meshEntry.clear() + else: + meshEntry = {} + xcord = [] + ycord = [] + meshEntryWildcard = '/##[ISA=ChemCompt]' + if modelRoot != '/': + meshEntryWildcard = modelRoot+meshEntryWildcard + for meshEnt in moose.wildcardFind(meshEntryWildcard): + mollist = [] + realist = [] + enzlist = [] + cplxlist = [] + tablist = [] + funclist = [] + + mol_cpl = moose.wildcardFind(meshEnt.path+'/##[ISA=PoolBase]') + funclist = moose.wildcardFind(meshEnt.path+'/##[ISA=Function]') + enzlist = moose.wildcardFind(meshEnt.path+'/##[ISA=EnzBase]') + realist = moose.wildcardFind(meshEnt.path+'/##[ISA=ReacBase]') + tablist = moose.wildcardFind(meshEnt.path+'/##[ISA=StimulusTable]') + if mol_cpl or funclist or enzlist or realist or tablist: + for m in mol_cpl: + if isinstance(moose.element(m.parent),moose.CplxEnzBase): + cplxlist.append(m) + objInfo = m.parent.path+'/info' + elif isinstance(moose.element(m),moose.PoolBase): + mollist.append(m) + objInfo =m.path+'/info' + if moose.exists(objInfo): + listOfitems[moose.element(moose.element(objInfo).parent)]={'x':xyPosition(objInfo,'x'),'y':xyPosition(objInfo,'y')} + + xcord.append(xyPosition(objInfo,'x')) + ycord.append(xyPosition(objInfo,'y')) + + getxyCord(xcord,ycord,funclist,listOfitems) + getxyCord(xcord,ycord,enzlist,listOfitems) + getxyCord(xcord,ycord,realist,listOfitems) + getxyCord(xcord,ycord,tablist,listOfitems) + + meshEntry[meshEnt] = {'enzyme':enzlist, + 'reaction':realist, + 'pool':mollist, + 'cplx':cplxlist, + 'table':tablist, + 'function':funclist + } + positionInfoExist = not(len(np.nonzero(xcord)[0]) == 0 \ + and len(np.nonzero(ycord)[0]) == 0) + if positionInfoExist: + xmin = min(xcord) + xmax = max(xcord) + ymin = min(ycord) + ymax = max(ycord) + return meshEntry,xmin,xmax,ymin,ymax,positionInfoExist,listOfitems + +def sizeHint(self): + return QtCore.QSize(800,400) + +def getxyCord(xcord,ycord,list1,listOfitems): + for item in list1: + # if isinstance(item,Function): + # objInfo = moose.element(item.parent).path+'/info' + # else: + # objInfo = item.path+'/info' + if not isinstance(item,moose.Function): + objInfo = item.path+'/info' + xcord.append(xyPosition(objInfo,'x')) + ycord.append(xyPosition(objInfo,'y')) + if moose.exists(objInfo): + listOfitems[moose.element(moose.element(objInfo).parent)]={'x':xyPosition(objInfo,'x'),'y':xyPosition(objInfo,'y')} + +def setupItem(modelPath,cntDict): + '''This function collects information of what is connected to what. \ + eg. substrate and product connectivity to reaction's and enzyme's \ + sumtotal connectivity to its pool are collected ''' + #print " setupItem" + sublist = [] + prdlist = [] + zombieType = ['ReacBase','EnzBase','Function','StimulusTable'] + for baseObj in zombieType: + path = '/##[ISA='+baseObj+']' + if modelPath != '/': + path = modelPath+path + if ( (baseObj == 'ReacBase') or (baseObj == 'EnzBase')): + for items in moose.wildcardFind(path): + sublist = [] + prdlist = [] + uniqItem,countuniqItem = countitems(items,'subOut') + subNo = uniqItem + for sub in uniqItem: + sublist.append((moose.element(sub),'s',countuniqItem[sub])) + + uniqItem,countuniqItem = countitems(items,'prd') + prdNo = uniqItem + if (len(subNo) == 0 or len(prdNo) == 0): + print ("Substrate Product is empty ",path, " ",items) + + for prd in uniqItem: + prdlist.append((moose.element(prd),'p',countuniqItem[prd])) + + if (baseObj == 'CplxEnzBase') : + uniqItem,countuniqItem = countitems(items,'toEnz') + for enzpar in uniqItem: + sublist.append((moose.element(enzpar),'t',countuniqItem[enzpar])) + + uniqItem,countuniqItem = countitems(items,'cplxDest') + for cplx in uniqItem: + prdlist.append((moose.element(cplx),'cplx',countuniqItem[cplx])) + + if (baseObj == 'EnzBase'): + uniqItem,countuniqItem = countitems(items,'enzDest') + for enzpar in uniqItem: + sublist.append((moose.element(enzpar),'t',countuniqItem[enzpar])) + cntDict[items] = sublist,prdlist + elif baseObj == 'Function': + for items in moose.wildcardFind(path): + sublist = [] + prdlist = [] + item = items.path+'/x[0]' + uniqItem,countuniqItem = countitems(item,'input') + for funcpar in uniqItem: + sublist.append((moose.element(funcpar),'sts',countuniqItem[funcpar])) + + uniqItem,countuniqItem = countitems(items,'valueOut') + for funcpar in uniqItem: + prdlist.append((moose.element(funcpar),'stp',countuniqItem[funcpar])) + cntDict[items] = sublist,prdlist + + # elif baseObj == 'Function': + # #ZombieSumFunc adding inputs + # inputlist = [] + # outputlist = [] + # funplist = [] + # nfunplist = [] + # for items in moose.wildcardFind(path): + # for funplist in moose.element(items).neighbors['valueOut']: + # for func in funplist: + # funcx = moose.element(items.path+'/x[0]') + # uniqItem,countuniqItem = countitems(funcx,'input') + # for inPut in uniqItem: + # inputlist.append((inPut,'st',countuniqItem[inPut])) + # cntDict[func] = inputlist + else: + for tab in moose.wildcardFind(path): + tablist = [] + uniqItem,countuniqItem = countitems(tab,'output') + for tabconnect in uniqItem: + tablist.append((moose.element(tabconnect),'tab',countuniqItem[tabconnect])) + cntDict[tab] = tablist + +def countitems(mitems,objtype): + items = [] + items = moose.element(mitems).neighbors[objtype] + uniqItems = set(items) + countuniqItems = Counter(items) + return(uniqItems,countuniqItems) \ No newline at end of file diff --git a/moose-core/python/moose/chemUtil/graphUtils.py b/moose-core/python/moose/chemUtil/graphUtils.py new file mode 100644 index 0000000000000000000000000000000000000000..1d9ee1d550d05576191ef141c963359dbeb52c2f --- /dev/null +++ b/moose-core/python/moose/chemUtil/graphUtils.py @@ -0,0 +1,77 @@ +import moose +import pygraphviz as pgv +#import networkx as nx + +def autoCoordinates(meshEntry,srcdesConnection): + positionInfo = {} + + if meshEntry: + #G = nx.Graph() + G = pgv.AGraph() + + for cmpt,memb in list(meshEntry.items()): + for enzObj in find_index(memb,'enzyme'): + #G.add_node(enzObj.path) + G.add_node(enzObj.path,label='',shape='ellipse',color='',style='filled',fontname='Helvetica',fontsize=12,fontcolor='blue') + for cmpt,memb in list(meshEntry.items()): + for poolObj in find_index(memb,'pool'): + #G.add_node(poolObj.path) + G.add_node(poolObj.path,label = poolObj.name,shape = 'box',color = '',style = 'filled',fontname = 'Helvetica',fontsize = 12,fontcolor = 'blue') + for cplxObj in find_index(memb,'cplx'): + G.add_node(cplxObj.path) + G.add_node(cplxObj.path,label = cplxObj.name,shape = 'box',color = '',style = 'filled',fontname = 'Helvetica',fontsize = 12,fontcolor = 'blue') + #G.add_edge((cplxObj.parent).path,cplxObj.path) + for reaObj in find_index(memb,'reaction'): + #G.add_node(reaObj.path) + G.add_node(reaObj.path,label='',shape='circle',color='') + + for inn,out in list(srcdesConnection.items()): + if (inn.className =='ZombieReac'): arrowcolor = 'green' + elif(inn.className =='ZombieEnz'): arrowcolor = 'red' + else: arrowcolor = 'blue' + if isinstance(out,tuple): + if len(out[0])== 0: + print( inn.className + ':' +inn.name + " doesn't have input message") + else: + for items in (items for items in out[0] ): + G.add_edge(moose.element(items[0]).path,inn.path) + if len(out[1]) == 0: + print(inn.className + ':' + inn.name + "doesn't have output mssg") + else: + for items in (items for items in out[1] ): + G.add_edge(inn.path,moose.element(items[0]).path) + elif isinstance(out,list): + if len(out) == 0: + print ("Func pool doesn't have sumtotal") + else: + for items in (items for items in out ): + G.add_edge(moose.element(items[0]).path,inn.path) + + # if int( nx.__version__.split( '.' )[-1] ) >= 11: + # position = nx.spring_layout( G ) + # else: + # position = nx.graphviz_layout(G, prog = 'dot') + # for item in position.items(): + # xy = item[1] + # ann = moose.Annotator(item[0]+'/info') + # ann.x = xy[0] + # ann.y = xy[1] + # positioninfo[(moose.element(item[0]))] ={'x':float(xy[0]),'y':float(xy[1])} + + G.layout() + for n in G.nodes(): + value = str(n.attr['pos']) + valuelist = (value.split(',')) + positionInfo[(moose.element(n))] ={'x':float(valuelist[0]),'y':float(valuelist[1])} + ann = moose.Annotator(moose.element(n).path+'/info') + ann.x = float(valuelist[0]) + ann.y = float(valuelist[1]) + + return(positionInfo) + +def find_index(value, key): + """ Value.get(key) to avoid expection which would raise if empty value in dictionary for a given key """ + if value.get(key) != None: + return value.get(key) + else: + raise ValueError('no dict with the key found') diff --git a/moose-core/python/moose/genesis/writeKkit.py b/moose-core/python/moose/genesis/writeKkit.py index 226fe2a319b000f199357f41a3fac3fdf0f85a73..ba71603148d5e5096a41e14499027a7a11b442c2 100644 --- a/moose-core/python/moose/genesis/writeKkit.py +++ b/moose-core/python/moose/genesis/writeKkit.py @@ -1,7 +1,7 @@ """ Chemical Signalling model loaded into moose can be save into Genesis-Kkit format """ __author__ = "Harsha Rani" -__copyright__ = "Copyright 2015, Harsha Rani and NCBS Bangalore" +__copyright__ = "Copyright 2017, Harsha Rani and NCBS Bangalore" __credits__ = ["NCBS Bangalore"] __license__ = "GNU GPL" __version__ = "1.0.0" @@ -11,12 +11,11 @@ __status__ = "Development" import sys import random -import moose -import numpy as np import re -from collections import Counter -import networkx as nx -import matplotlib +import matplotlib +import moose +from moose.chemUtil.chemConnectUtil import * +from moose.chemUtil.graphUtils import * GENESIS_COLOR_SEQUENCE = ((248, 0, 255), (240, 0, 255), (232, 0, 255), (224, 0, 255), (216, 0, 255), (208, 0, 255), (200, 0, 255), (192, 0, 255), (184, 0, 255), (176, 0, 255), (168, 0, 255), (160, 0, 255), (152, 0, 255), (144, 0, 255), @@ -38,7 +37,7 @@ GENESIS_COLOR_SEQUENCE = ((248, 0, 255), (240, 0, 255), (232, 0, 255), (224, 0, #Todo : To be written # --StimulusTable -def mooseWriteKkit( modelpath, filename,sceneitems=None): +def mooseWriteKkit( modelpath, filename, sceneitems={}): if filename.rfind('.') != -1: filename = filename[:filename.rfind('.')] else: @@ -47,45 +46,37 @@ def mooseWriteKkit( modelpath, filename,sceneitems=None): global NA NA = 6.0221415e23 global cmin,cmax,xmin,xmax,ymin,ymax - cmin = 0 - cmax = 1 - xmin = 0 - xmax = 1 - ymin = 0 - ymax = 1 - + cmin, xmin, ymin = 0, 0, 0 + cmax, xmax, ymax = 1, 1, 1 + compt = moose.wildcardFind(modelpath+'/##[ISA=ChemCompt]') maxVol = estimateDefaultVol(compt) - f = open(filename, 'w') - - if sceneitems == None: - srcdesConnection = {} - setupItem(modelpath,srcdesConnection) - meshEntry = setupMeshObj(modelpath) - cmin,cmax,sceneitems = autoCoordinates(meshEntry,srcdesConnection) - for k,v in list(sceneitems.items()): - v = sceneitems[k] - x1 = calPrime(v['x']) - y1 = calPrime(v['y']) - sceneitems[k]['x'] = x1 - sceneitems[k]['y'] = y1 - else: - cs, xcord, ycord = [], [] ,[] - for k,v in list(sceneitems.items()): - xcord.append(v['x']) - cs.append(v['x']) - ycord.append(v['y']) - cs.append(v['y']) - xmin = min(xcord) - xmax = max(xcord) - ymin = min(ycord) - ymax = max(ycord) - - cmin = min(cs) - cmax = max(cs) - writeHeader (f,maxVol) - - if (compt > 0): + positionInfoExist = True + if compt: + if bool(sceneitems): + cmin,cmax,xmin1,xmax1,ymin1,ymax1 = findMinMax(sceneitems) + elif not bool(sceneitems): + srcdesConnection = {} + setupItem(modelpath,srcdesConnection) + meshEntry,xmin,xmax,ymin,ymax,positionInfoExist,sceneitems = setupMeshObj(modelpath) + if not positionInfoExist: + #cmin,cmax,sceneitems = autoCoordinates(meshEntry,srcdesConnection) + sceneitems = autoCoordinates(meshEntry,srcdesConnection) + + if not positionInfoExist: + # if position are not from kkit, then zoom factor is applied while + # writing to genesis. Like if position is from pyqtSceneItem or auto-coordinates + cmin,cmax,xmin1,xmax1,ymin1,ymax1 = findMinMax(sceneitems) + for k,v in list(sceneitems.items()): + anno = moose.element(k.path+'/info') + x1 = calPrime(v['x']) + y1 = calPrime(v['y']) + sceneitems[k]['x'] = x1 + sceneitems[k]['y'] = y1 + + f = open(filename, 'w') + writeHeader (f,maxVol) + gtId_vol = writeCompartment(modelpath,compt,f) writePool(modelpath,f,gtId_vol,sceneitems) reacList = writeReac(modelpath,f,sceneitems) @@ -114,188 +105,36 @@ def mooseWriteKkit( modelpath, filename,sceneitems=None): writeFooter1(f) writeNotes(modelpath,f) writeFooter2(f) + print('Written to file '+filename) return True else: print(("Warning: writeKkit:: No model found on " , modelpath)) return False -def calPrime(x): - prime = int((20*(float(x-cmin)/float(cmax-cmin)))-10) - return prime - -def setupItem(modelPath,cntDict): - '''This function collects information of what is connected to what. \ - eg. substrate and product connectivity to reaction's and enzyme's \ - sumtotal connectivity to its pool are collected ''' - #print " setupItem" - sublist = [] - prdlist = [] - zombieType = ['ReacBase','EnzBase','Function','StimulusTable'] - for baseObj in zombieType: - path = '/##[ISA='+baseObj+']' - if modelPath != '/': - path = modelPath+path - if ( (baseObj == 'ReacBase') or (baseObj == 'EnzBase')): - for items in moose.wildcardFind(path): - sublist = [] - prdlist = [] - uniqItem,countuniqItem = countitems(items,'subOut') - subNo = uniqItem - for sub in uniqItem: - sublist.append((moose.element(sub),'s',countuniqItem[sub])) - - uniqItem,countuniqItem = countitems(items,'prd') - prdNo = uniqItem - if (len(subNo) == 0 or len(prdNo) == 0): - print("Substrate Product is empty ",path, " ",items) - - for prd in uniqItem: - prdlist.append((moose.element(prd),'p',countuniqItem[prd])) - - if (baseObj == 'CplxEnzBase') : - uniqItem,countuniqItem = countitems(items,'toEnz') - for enzpar in uniqItem: - sublist.append((moose.element(enzpar),'t',countuniqItem[enzpar])) - - uniqItem,countuniqItem = countitems(items,'cplxDest') - for cplx in uniqItem: - prdlist.append((moose.element(cplx),'cplx',countuniqItem[cplx])) - - if (baseObj == 'EnzBase'): - uniqItem,countuniqItem = countitems(items,'enzDest') - for enzpar in uniqItem: - sublist.append((moose.element(enzpar),'t',countuniqItem[enzpar])) - cntDict[items] = sublist,prdlist - elif baseObj == 'Function': - for items in moose.wildcardFind(path): - sublist = [] - prdlist = [] - item = items.path+'/x[0]' - uniqItem,countuniqItem = countitems(item,'input') - for funcpar in uniqItem: - sublist.append((moose.element(funcpar),'sts',countuniqItem[funcpar])) - - uniqItem,countuniqItem = countitems(items,'valueOut') - for funcpar in uniqItem: - prdlist.append((moose.element(funcpar),'stp',countuniqItem[funcpar])) - cntDict[items] = sublist,prdlist - else: - for tab in moose.wildcardFind(path): - tablist = [] - uniqItem,countuniqItem = countitems(tab,'output') - for tabconnect in uniqItem: - tablist.append((moose.element(tabconnect),'tab',countuniqItem[tabconnect])) - cntDict[tab] = tablist -def countitems(mitems,objtype): - items = [] - items = moose.element(mitems).neighbors[objtype] - uniqItems = set(items) - countuniqItems = Counter(items) - return(uniqItems,countuniqItems) - -def setupMeshObj(modelRoot): - ''' Setup compartment and its members pool,reaction,enz cplx under self.meshEntry dictionaries \ - self.meshEntry with "key" as compartment, - value is key2:list where key2 represents moose object type,list of objects of a perticular type - e.g self.meshEntry[meshEnt] = { 'reaction': reaction_list,'enzyme':enzyme_list,'pool':poollist,'cplx': cplxlist } - ''' - meshEntry = {} - if meshEntry: - meshEntry.clear() - else: - meshEntry = {} - meshEntryWildcard = '/##[ISA=ChemCompt]' - if modelRoot != '/': - meshEntryWildcard = modelRoot+meshEntryWildcard - for meshEnt in moose.wildcardFind(meshEntryWildcard): - mollist = [] - cplxlist = [] - mol_cpl = moose.wildcardFind(meshEnt.path+'/##[ISA=PoolBase]') - funclist = moose.wildcardFind(meshEnt.path+'/##[ISA=Function]') - enzlist = moose.wildcardFind(meshEnt.path+'/##[ISA=EnzBase]') - realist = moose.wildcardFind(meshEnt.path+'/##[ISA=ReacBase]') - tablist = moose.wildcardFind(meshEnt.path+'/##[ISA=StimulusTable]') - if mol_cpl or funclist or enzlist or realist or tablist: - for m in mol_cpl: - if isinstance(moose.element(m.parent),moose.CplxEnzBase): - cplxlist.append(m) - elif isinstance(moose.element(m),moose.PoolBase): - mollist.append(m) - - meshEntry[meshEnt] = {'enzyme':enzlist, - 'reaction':realist, - 'pool':mollist, - 'cplx':cplxlist, - 'table':tablist, - 'function':funclist - } - return(meshEntry) -def autoCoordinates(meshEntry,srcdesConnection): - G = nx.Graph() - for cmpt,memb in list(meshEntry.items()): - for enzObj in find_index(memb,'enzyme'): - #G.add_node(enzObj.path) - G.add_node(enzObj.path,label=enzObj.name,shape='ellipse',color='',style='filled',fontname='Helvetica',fontsize=12,fontcolor='blue') - for cmpt,memb in list(meshEntry.items()): - for poolObj in find_index(memb,'pool'): - #G.add_node(poolObj.path) - G.add_node(poolObj.path,label = poolObj.name,shape = 'box',color = '',style = 'filled',fontname = 'Helvetica',fontsize = 12,fontcolor = 'blue') - for cplxObj in find_index(memb,'cplx'): - pass - #G.add_node(cplxObj.path) - #G.add_edge((cplxObj.parent).path,cplxObj.path) - for reaObj in find_index(memb,'reaction'): - #G.add_node(reaObj.path) - G.add_node(reaObj.path,label=reaObj.name,shape='',color='') - for funcObj in find_index(memb,'function'): - G.add_node(poolObj.path,label = funcObj.name,shape = 'box',color = 'red',style = 'filled',fontname = 'Helvetica',fontsize = 12,fontcolor = 'blue') - - - for inn,out in list(srcdesConnection.items()): - if (inn.className =='ZombieReac'): arrowcolor = 'green' - elif(inn.className =='ZombieEnz'): arrowcolor = 'red' - else: arrowcolor = 'blue' - if isinstance(out,tuple): - if len(out[0])== 0: - print(inn.className + ':' +inn.name + " doesn't have input message") - else: - for items in (items for items in out[0] ): - G.add_edge(moose.element(items[0]).path,inn.path) - if len(out[1]) == 0: - print(inn.className + ':' + inn.name + "doesn't have output mssg") - else: - for items in (items for items in out[1] ): - G.add_edge(inn.path,moose.element(items[0]).path) - elif isinstance(out,list): - if len(out) == 0: - print("Func pool doesn't have sumtotal") - else: - for items in (items for items in out ): - G.add_edge(moose.element(items[0]).path,inn.path) - - position = nx.graphviz_layout(G, prog = 'dot') - if int( nx.__version__.split( '.' )[-1] ) >= 11: - position = nx.spring_layout( G ) - - #agraph = nx.to_agraph(G) - #agraph.draw("writetogenesis.png", format = 'png', prog = 'dot') - sceneitems = {} +def findMinMax(sceneitems): + cmin = 0.0 + cmax = 1.0 + xmin,xymin = 0.0,0.0 + xmax,xymax = 1.0,1.0 xycord = [] - - for key,value in list(position.items()): - xycord.append(value[0]) - xycord.append(value[1]) - sceneitems[moose.element(key)] = {'x':value[0],'y':value[1]} + xcord = [] + ycord = [] + for k,v in list(sceneitems.items()): + xycord.append(v['x']) + xycord.append(v['y']) + xcord.append(v['x']) + ycord.append(v['y']) + xmin = min(xcord) + xmax = max(xcord) + ymin = min(ycord) + ymax = max(ycord) cmin = min(xycord) cmax = max(xycord) - return cmin,cmax,sceneitems + return cmin,cmax,xmin,xmax,ymin,ymax -def find_index(value, key): - """ Value.get(key) to avoid expection which would raise if empty value in dictionary for a given key """ - if value.get(key) != None: - return value.get(key) - else: - raise ValueError('no dict with the key found') +def calPrime(x): + prime = int((20*(float(x-cmin)/float(cmax-cmin)))-10) + return prime def storeCplxEnzMsgs( enz, f ): for sub in enz.neighbors["subOut"]: @@ -409,9 +248,7 @@ def nearestColorIndex(color, color_sequence): #Trying to find the index to closest color map from the rainbow pickle file for matching the Genesis color map distance = [ (color[0] - temp[0]) ** 2 + (color[1] - temp[1]) ** 2 + (color[2] - temp[2]) ** 2 for temp in color_sequence] - minindex = 0 - for i in range(1, len(distance)): if distance[minindex] > distance[i] : minindex = i @@ -500,11 +337,23 @@ def storePlotMsgs( tgraphs,f): for graph in tgraphs: slash = graph.path.find('graphs') if not slash > -1: - slash = graph.path.find('graph_0') + slash = graph.path.find('graph') if slash > -1: - conc = graph.path.find('conc') - if conc > -1 : - tabPath = graph.path[slash:len(graph.path)] + foundConc = True + if not ( (graph.path.find('conc1') > -1 ) or + (graph.path.find('conc2') > -1 ) or + (graph.path.find('conc3') > -1 ) or + (graph.path.find('conc4') > -1) ): + foundConc = False + + #conc = graph.path.find('conc') + # if conc > -1 : + # tabPath = graph.path[slash:len(graph.path)] + # else: + # slash1 = graph.path.find('/',slash) + # tabPath = "/graphs/conc1" +graph.path[slash1:len(graph.path)] + if foundConc == True: + tabPath = "/"+graph.path[slash:len(graph.path)] else: slash1 = graph.path.find('/',slash) tabPath = "/graphs/conc1" +graph.path[slash1:len(graph.path)] @@ -527,14 +376,21 @@ def writeplot( tgraphs,f ): for graphs in tgraphs: slash = graphs.path.find('graphs') if not slash > -1: - slash = graphs.path.find('graph_0') + slash = graphs.path.find('graph') if slash > -1: - conc = graphs.path.find('conc') - if conc > -1 : + foundConc = True + if not ( (graphs.path.find('conc1') > -1 ) or + (graphs.path.find('conc2') > -1 ) or + (graphs.path.find('conc3') > -1 ) or + (graphs.path.find('conc4') > -1) ): + foundConc = False + if foundConc == True: tabPath = "/"+graphs.path[slash:len(graphs.path)] else: slash1 = graphs.path.find('/',slash) tabPath = "/graphs/conc1" +graphs.path[slash1:len(graphs.path)] + + if len(moose.element(graphs).msgOut): poolPath = (moose.element(graphs).msgOut)[0].e2.path poolEle = moose.element(poolPath) @@ -586,7 +442,7 @@ def writePool(modelpath,f,volIndex,sceneitems): color = getRandColor() if textcolor == "" or textcolor == " ": textcolor = getRandColor() - #print " trimPath",trimPath(p) + #print " trimPath",trimPath(p) f.write("simundump kpool /kinetics/" + trimPath(p) + " 0 " + str(p.diffConst) + " " + str(0) + " " + @@ -765,7 +621,7 @@ if __name__ == "__main__": filename = sys.argv[1] modelpath = filename[0:filename.find('.')] moose.loadModel(filename,'/'+modelpath,"gsl") - output = modelpath+"_4mmoose.g" + output = modelpath+"_.g" written = write('/'+modelpath,output) if written: print((" file written to ",output)) diff --git a/moose-core/python/moose/merge.py b/moose-core/python/moose/merge.py deleted file mode 100644 index 70f37b24530de75e7aaf80dd7c30b152e55198b5..0000000000000000000000000000000000000000 --- a/moose-core/python/moose/merge.py +++ /dev/null @@ -1,248 +0,0 @@ -from __future__ import print_function - -#******************************************************************* -# * File: merge.py -# * Description: -# * Author: HarshaRani -# * E-mail: hrani@ncbs.res.in -# ********************************************************************/ -# ********************************************************************** -#** This program is part of 'MOOSE', the -#** Messaging Object Oriented Simulation Environment, -#** also known as GENESIS 3 base code. -#** copyright (C) 2003-2016 Upinder S. Bhalla. and NCBS -#Created : Friday June 13 12:19:00 2016(+0530) -#Version -#Last-Updated:Tuesday June 21 17.48.37 -# By: Harsha -#**********************************************************************/ - -# This program is used to merge models -# -- Model B is merged to A -# -- Rules are -# ## check for compartment in B exist in A -# *** If compartment from model B doesn't exist in model A, then copy entire compartment from B to A -# *** If compartment from model B exist in model A -# ^^^^ check the volume of compartment of B as compared to A -# !!!!! If same, then copy all the moose object which doesn't exist in A and reference to both mooseId -# which is used for connecting objects. -# -# ^^^^ If volume of compartment of model B is different then change the volume of compartment of model B same as A -# !!!! If same, then copy all the moose object which doesn't exist in A and reference to both mooseId -# which is used for connecting objects - -# &&&&&& copying pools from compartment from B to A is easy -# &&&&& but while copying reaction and enzyme check if the same reaction exist -# -- if same reaction name exists -# -- Yes -# 1) check if substrate and product are same as compartment from model B to compartment modelA -# --yes do nothing -# --No then duplicate the reaction in model A and connect the message of approraite sub/prd -# and warn the user the same reaction name with different sub or product -# -- No -# copy the reaction - -import sys -from . import _moose as moose - -def merge(A,B): - #load models into moose and solver's are deleted - apath = loadModels(A) - bpath = loadModels(B) - - comptAdict = comptList(apath) - comptBdict = comptList(bpath) - - for key in list(comptBdict.keys()): - - if key not in comptAdict: - # comptBdict[key] - compartment from model B which does not exist in model A - moose.copy(comptBdict[key],moose.element(apath)) - else: - war_msg = "" - copied , duplicated = [],[] - - if not (comptAdict[key].volume == comptBdict[key].volume): - #change volume of ModelB same as ModelA - volumeA = comptAdict[key].volume - comptBdict[key].volume = volumeA - - #Merging Pool - poolName_a = [] - poolListina = moose.wildcardFind(comptAdict[key].path+'/##[ISA=PoolBase]') - if poolListina: - neutral_compartment = findgroup_compartment(poolListina[0]) - for pl in poolListina: - poolName_a.append(pl.name) - - for pool in moose.wildcardFind(comptBdict[key].path+'/##[ISA=PoolBase]'): - if not pool.name in poolName_a : - #pool has compartment(assume its direct pool ),so copy under kinetics - if ((pool.parent).className == "CubeMesh" or (pool.parent).className == "Neutral"): - c = moose.copy(pool,neutral_compartment) - copied.append(c) - elif (pool.parent).className == "ZombieEnz" or (pool.parent).className == "Enz": - print(" Pool not in modelA but enz parent",pool.name) - pass - else: - print(" Check this pool, parent which doesn't have ChemCompt or Enz", end=' ') - - #Mergering StimulusTable - stimName_a = [] - stimListina = moose.wildcardFind(comptAdict[key].path+'/##[ISA=StimulusTable]') - if stimListina: - neutral_compartment = findgroup_compartment(stimListina[0]) - for st in stimListina: - stimName_a.append(st.name) - for stb in moose.wildcardFind(comptBdict[key].path+'/##[ISA=StimulusTable]'): - if stb.name in stimName_a: - sA = comptAdict[key].path+'/'+stb.name - sB = comptBdict[key].path+'/'+stb.name - stAOutput = subprdList(sA,"output") - stBOutput = subprdList(sB,"output") - sas = set(stAOutput) - sbs = set(stBOutput) - uniq = (sas.union(sbs) - sas.intersection(sbs)) - for u in uniq: - if u not in stAOutput: - src = moose.element(sA) - des = moose.element(comptAdict[key].path+'/'+u) - moose.connect(src,'output',des,'setConcInit') - else: - st1 = moose.StimulusTable(comptAdict[key].path+'/'+stb.name) - for sb in sbs: - des = moose.element(comptAdict[key].path+'/'+sb) - moose.connect(st1,'output',des,'setConcInit') - #Mergering Reaction - reacName_a = [] - reacListina = moose.wildcardFind(comptAdict[key].path+'/##[ISA=ReacBase]') - if reacListina: - neutral_compartment = findgroup_compartment(poolListina[0]) - - for re in reacListina: - reacName_a.append(re.name) - for reac in moose.wildcardFind(comptBdict[key].path+'/##[ISA=ReacBase]'): - if reac.name in reacName_a: - - rA = comptAdict[key].path+'/'+reac.name - rB = comptBdict[key].path+'/'+reac.name - - rAsubname,rBsubname,rAprdname,rBprdname = [], [], [], [] - rAsubname = subprdList(rA,"sub") - rBsubname = subprdList(rB,"sub") - rAprdname = subprdList(rA,"prd") - rBprdname = subprdList(rB,"prd") - - aS = set(rAsubname) - bS = set(rBsubname) - aP = set(rAprdname) - bP = set(rBprdname) - - hasSameSub,hasSamePrd = True,True - hassameSlen,hassamePlen = False,False - - hasSameSub = not (len (aS.union(bS) - aS.intersection(bS))) - hasSamePrd = not (len(aP.union(bP) - aP.intersection(bP))) - - hassameSlen = ( len(rAsubname) == len(rBsubname)) - hassamePlen = ( len(rAprdname) == len(rBprdname)) - - - if not all((hasSameSub,hasSamePrd,hassameSlen,hassamePlen)): - war_msg = war_msg +"Reaction \""+reac.name+ "\" also contains in modelA but with different" - if not all((hasSameSub,hassameSlen)): - war_msg = war_msg+ " substrate " - - if not all((hasSamePrd, hassamePlen)): - war_msg = war_msg+ " product" - war_msg = war_msg +", reaction is duplicated in modelA. \nModeler should decide to keep or remove this reaction" - - reac = moose.Reac(comptAdict[key].path+'/'+reac.name+"_duplicated") - mooseConnect(comptAdict[key].path,reac,rBsubname,"sub") - mooseConnect(comptAdict[key].path,reac,rBprdname,"prd") - - duplicated.append(reac) - - elif not reac.name in reacName_a : - #reac has compartment(assume its direct reac ),so copy under kinetics - if ((reac.parent).className == "CubeMesh" or (reac.parent).className == "Neutral"): - c = moose.copy(reac,neutral_compartment) - copied.append(c) - - rBsubname, rBprdname = [],[] - rBsubname = subprdList(reac,"sub") - rBprdname = subprdList(reac,"prd") - mooseConnect(comptAdict[key].path,reac,rBsubname,"sub") - mooseConnect(comptAdict[key].path,reac,rBprdname,"prd") - print("\ncopied: ", copied, \ - "\n\nDuplicated: ",duplicated, \ - "\n\nwarning: ",war_msg) - return copied,duplicated,war_msg -def loadModels(filename): - apath = '/' - - apath = filename[filename.rfind('/'): filename.rfind('.')] - - moose.loadModel(filename,apath) - - #Solvers are removed - deleteSolver(apath) - return apath - -def comptList(modelpath): - comptdict = {} - for ca in moose.wildcardFind(modelpath+'/##[ISA=ChemCompt]'): - comptdict[ca.name] = ca - return comptdict - -def mooseConnect(modelpath,reac,spList,sptype): - for spl in spList: - spl_id = moose.element(modelpath+'/'+spl) - reac_id = moose.element(modelpath+'/'+reac.name) - if sptype == "sub": - m = moose.connect( reac_id, "sub", spl_id, 'reac' ) - else: - moose.connect( reac_id, "prd", spl_id, 'reac' ) - -def deleteSolver(modelRoot): - compts = moose.wildcardFind(modelRoot+'/##[ISA=ChemCompt]') - for compt in compts: - if moose.exists(compt.path+'/stoich'): - st = moose.element(compt.path+'/stoich') - st_ksolve = st.ksolve - moose.delete(st) - if moose.exists((st_ksolve).path): - moose.delete(st_ksolve) - -def subprdList(reac,subprd): - rtype = moose.element(reac).neighbors[subprd] - rname = [] - for rs in rtype: - rname.append(rs.name) - return rname - -def findCompartment(element): - while not mooseIsInstance(element,["CubeMesh","CyclMesh"]): - element = element.parent - return element -def mooseIsInstance(element, classNames): - #print classNames - #print moose.element(element).__class__.__name__ in classNames - return moose.element(element).__class__.__name__ in classNames - -def findgroup_compartment(element): - #Try finding Group which is Neutral but it should exist before Compartment, if Compartment exist then stop at Comparment - while not mooseIsInstance(element,"Neutral"): - if mooseIsInstance(element,["CubeMesh","CyclMesh"]): - return element - element = element.parent - return element -if __name__ == "__main__": - - model1 = '/home/harsha/genesis_files/gfile/acc12.g' - model2 = '/home/harsha/Trash/acc12_withadditionPool.g' - #model1 = '/home/harsha/Trash/modelA.g' - #model2 = '/home/harsha/Trash/modelB.g' - model1 = '/home/harsha/genesis_files/gfile/acc44.g' - model2 = '/home/harsha/genesis_files/gfile/acc45.g' - mergered = merge(model1,model2) diff --git a/moose-core/python/moose/merge/merge.py b/moose-core/python/moose/merge/merge.py new file mode 100644 index 0000000000000000000000000000000000000000..f76ca494b3daffccc1fa6318dd9783564c3a09d1 --- /dev/null +++ b/moose-core/python/moose/merge/merge.py @@ -0,0 +1,609 @@ +#******************************************************************* +# * File: merge.py +# * Description: +# * Author: HarshaRani +# * E-mail: hrani@ncbs.res.in +# ********************************************************************/ +# ********************************************************************** +#** This program is part of 'MOOSE', the +#** Messaging Object Oriented Simulation Environment, +#** also known as GENESIS 3 base code. +#** copyright (C) 2003-2017 Upinder S. Bhalla. and NCBS +#Created : Friday Dec 16 23:19:00 2016(+0530) +#Version +#Last-Updated: Thursday Jan 12 17:30:33 2017(+0530) +# By: Harsha +#**********************************************************************/ + +# This program is used to merge models +# -- Model B is merged to modelA +#Rules are +#--In this models are mergered at group level (if exists) + + + +import sys +import os +#from . import _moose as moose +import moose +import mtypes + +from moose.chemUtil.chemConnectUtil import * +from moose.chemUtil.graphUtils import * + +def mergeChemModel(A,B): + """ Merges two model or the path """ + + modelA,loadedA = loadModels(A) + modelB,loadedB = loadModels(B) + if not loadedA or not loadedB: + if not loadedA: + modelB = moose.Shell('/') + if not loadedB: + modelA = moose.Shell('/') + else: + directory, bfname = os.path.split(B) + global grpNotcopiedyet,poolListina + poolListina = {} + grpNotcopiedyet = [] + dictComptA = dict( [ (i.name,i) for i in moose.wildcardFind(modelA+'/##[ISA=ChemCompt]') ] ) + dictComptB = dict( [ (i.name,i) for i in moose.wildcardFind(modelB+'/##[ISA=ChemCompt]') ] ) + + poolNotcopiedyet = [] + + for key in list(dictComptB.keys()): + if key not in dictComptA: + # if compartmentname from modelB does not exist in modelA, then copy + copy = moose.copy(dictComptB[key],moose.element(modelA)) + else: + #if compartmentname from modelB exist in modelA, + #volume is not same, then change volume of ModelB same as ModelA + if abs(dictComptA[key].volume - dictComptB[key].volume): + #hack for now + while (abs(dictComptA[key].volume - dictComptB[key].volume) != 0.0): + dictComptB[key].volume = float(dictComptA[key].volume) + dictComptA = dict( [ (i.name,i) for i in moose.wildcardFind(modelA+'/##[ISA=ChemCompt]') ] ) + + #Mergering pool + poolMerge(dictComptA[key],dictComptB[key],poolNotcopiedyet) + + if grpNotcopiedyet: + # objA = moose.element(comptA).parent.name + # if not moose.exists(objA+'/'+comptB.name+'/'+bpath.name): + # print bpath + # moose.copy(bpath,moose.element(objA+'/'+comptB.name)) + pass + + comptAdict = comptList(modelA) + poolListina = {} + poolListina = updatePoolList(comptAdict) + funcNotallowed = [] + R_Duplicated, R_Notcopiedyet,R_Daggling = [], [], [] + E_Duplicated, E_Notcopiedyet,E_Daggling = [], [], [] + for key in list(dictComptB.keys()): + funcNotallowed = [] + funcNotallowed = functionMerge(dictComptA,dictComptB,key) + + poolListina = updatePoolList(dictComptA) + R_Duplicated,R_Notcopiedyet,R_Daggling = reacMerge(dictComptA,dictComptB,key,poolListina) + + poolListina = updatePoolList(dictComptA) + E_Duplicated,E_Notcopiedyet,E_Daggling = enzymeMerge(dictComptA,dictComptB,key,poolListina) + + print("\n Model is merged to %s" %modelA) + + if funcNotallowed: + print( "\nPool already connected to a function, this function is not to connect to same pool, since no two function are allowed to connect to same pool:") + for fl in list(funcNotallowed): + print("\t [Pool]: %s [Function]: %s \n" %(str(fl.parent.name), str(fl.path))) + + if R_Duplicated or E_Duplicated: + print ("Reaction / Enzyme are Duplicate" + "\n 1. The once whoes substrate / product names are different for a give reaction name " + "\n 2. its compartment to which it belongs to may be is different" + "\n Models have to decide to keep or delete these reaction/enzyme") + if E_Duplicated: + print("Reaction: ") + for rd in list(R_Duplicated): + print ("%s " %str(rd.name)) + + if E_Duplicated: + print ("Enzyme:") + for ed in list(E_Duplicated): + print ("%s " %str(ed.name)) + + if R_Notcopiedyet or E_Notcopiedyet: + + print ("\nThese Reaction/Enzyme in model are not dagging but while copying the associated substrate or product is missing") + if R_Notcopiedyet: + print("Reaction: ") + for rd in list(R_Notcopiedyet): + print ("%s " %str(rd.name)) + if E_Notcopiedyet: + print ("Enzyme:") + for ed in list(E_Notcopiedyet): + print ("%s " %str(ed.name)) + + if R_Daggling or E_Daggling: + print ("\n Daggling reaction/enzyme are not not allowed in moose, these are not merged") + if R_Daggling: + print("Reaction: ") + for rd in list(R_Daggling): + print ("%s " %str(rd.name)) + if E_Daggling: + print ("Enzyme:") + for ed in list(E_Daggling): + print ("%s " %str(ed.name)) + +def functionMerge(comptA,comptB,key): + funcNotallowed = [] + comptApath = moose.element(comptA[key]).path + comptBpath = moose.element(comptB[key]).path + funcListina = moose.wildcardFind(comptApath+'/##[ISA=PoolBase]') + funcListinb = moose.wildcardFind(comptBpath+'/##[ISA=Function]') + objA = moose.element(comptApath).parent.name + objB = moose.element(comptBpath).parent.name + #For function mergering pool name is taken as reference + funcNotcopiedyet = [] + for fb in funcListinb: + + if fb.parent.className in ['ZombiePool','Pool','ZombieBufPool','BufPool']: + objA = moose.element(comptApath).parent.name + fbpath = fb.path + #funcpath = fbpath[fbpath.find(findCompartment(fb).name)-1:len(fbpath)] + funcparentB = fb.parent.path + funcpath = fbpath.replace(objB,objA) + funcparentA = funcparentB.replace(objB,objA) + tt = moose.element(funcparentA).neighbors['setN'] + if tt: + funcNotallowed.append(fb) + else: + if len(moose.element(fb.path+'/x').neighbors["input"]): + #inputB = moose.element(fb.path+'/x').neighbors["input"] + inputB = subprdList(moose.element(fb.path+'/x'),"input") + inputB_expr = fb.expr + if moose.exists(funcpath): + #inputA = moose.element(objA+funcpath+'/x').neighbors["input"] + inputA = subprdList(moose.element(funcpath+'/x'),"input") + inputA_expr = moose.element(funcpath).expr + hassameExpr = False + if inputA_expr == inputB_expr: + hassameExpr = True + hassameLen,hassameS,hassameVols = same_len_name_vol(inputA,inputB) + if not all((hassameLen,hassameS,hassameVols,hassameExpr)): + fb.name = fb.name+'_duplicatedF' + createFunction(fb,inputB,objB,objA) + else: + #function doesnot exist then copy + if len(inputB): + volinput = [] + for inb in inputB: + volinput.append(findCompartment(moose.element(inb)).volume) + if len(set(volinput)) == 1: + # If all the input connected belongs to one compartment then copy + createFunction(fb,inputB,objB,objA) + else: + # moose doesn't allow function input to come from different compartment + funcNotallowed.append(fb) + return funcNotallowed + +def createFunction(fb,inputB,objB,objA): + fapath1 = fb.path.replace(objB,objA) + fapath = fapath1.replace('[0]','') + if not moose.exists(fapath): + # if fb.parent.className in ['CubeMesh','CyclMesh']: + # des = moose.Function('/'+objA+'/'+fb.parent.name+'/'+fb.name) + # elif fb.parent.className in ['Pool','ZombiePool','BufPool','ZombieBufPool']: + # for akey in list(poolListina[findCompartment(fb).name]): + # if fb.parent.name == akey.name: + # des = moose.Function(akey.path+'/'+fb.name) + des = moose.Function(fapath) + moose.connect(des, 'valueOut', moose.element(fapath).parent,'setN' ) + for src in inputB: + pool = ((src.path).replace(objB,objA)).replace('[0]','') + numVariables = des.numVars + expr = "" + expr = (des.expr+'+'+'x'+str(numVariables)) + expr = expr.lstrip("0 +") + expr = expr.replace(" ","") + des.expr = expr + moose.connect( pool, 'nOut', des.x[numVariables], 'input' ) + + #if fb.expr != des.expr: + # print "Function ",des, " which is duplicated from modelB, expression is different, this is tricky in moose to know what those constants are connected to " + # print "ModelB ", fb, fb.expr, "\nModelA ",des, des.expr + +def comptList(modelpath): + comptdict = {} + for ca in moose.wildcardFind(modelpath+'/##[ISA=ChemCompt]'): + comptdict[ca.name] = ca + return comptdict + +def loadModels(filename): + """ load models into moose if file, if moosepath itself it passes back the path and + delete solver if exist """ + + modelpath = '/' + loaded = False + + if os.path.isfile(filename) : + modelpath = filename[filename.rfind('/'): filename.rfind('.')] + ext = os.path.splitext(filename)[1] + filename = filename.strip() + modeltype = mtypes.getType(filename) + subtype = mtypes.getSubtype(filename, modeltype) + if subtype == 'kkit' or modeltype == "cspace": + moose.loadModel(filename,modelpath) + loaded = True + + elif subtype == 'sbml': + #moose.ReadSBML() + pass + else: + print("This file is not supported for mergering") + modelpath = moose.Shell('/') + elif moose.exists(filename): + modelpath = filename + loaded = True + ## default is 'ee' solver while loading the model using moose.loadModel, + ## yet deleteSolver is called just to be assured + if loaded: + deleteSolver(modelpath) + + return modelpath,loaded + +def deleteSolver(modelRoot): + compts = moose.wildcardFind(modelRoot+'/##[ISA=ChemCompt]') + for compt in compts: + if moose.exists(compt.path+'/stoich'): + st = moose.element(compt.path+'/stoich') + st_ksolve = st.ksolve + moose.delete(st) + if moose.exists((st_ksolve).path): + moose.delete(st_ksolve) + +def poolMerge(comptA,comptB,poolNotcopiedyet): + + aCmptGrp = moose.wildcardFind(comptA.path+'/#[TYPE=Neutral]') + aCmptGrp = aCmptGrp +(moose.element(comptA.path),) + + bCmptGrp = moose.wildcardFind(comptB.path+'/#[TYPE=Neutral]') + bCmptGrp = bCmptGrp +(moose.element(comptB.path),) + + objA = moose.element(comptA.path).parent.name + objB = moose.element(comptB.path).parent.name + for bpath in bCmptGrp: + grp_cmpt = ((bpath.path).replace(objB,objA)).replace('[0]','') + if moose.exists(grp_cmpt) : + if moose.element(grp_cmpt).className != bpath.className: + grp_cmpt = grp_cmpt+'_grp' + bpath.name = bpath.name+"_grp" + l = moose.Neutral(grp_cmpt) + else: + moose.Neutral(grp_cmpt) + + apath = moose.element(bpath.path.replace(objB,objA)) + + bpoollist = moose.wildcardFind(bpath.path+'/#[ISA=PoolBase]') + apoollist = moose.wildcardFind(apath.path+'/#[ISA=PoolBase]') + for bpool in bpoollist: + if bpool.name not in [apool.name for apool in apoollist]: + copied = copy_deleteUnlyingPoolObj(bpool,apath) + if copied == False: + #hold it for later, this pool may be under enzyme, as cplx + poolNotcopiedyet.append(bpool) + +def copy_deleteUnlyingPoolObj(pool,path): + # check if this pool is under compartement or under enzyme?(which is enzyme_cplx) + # if enzyme_cplx then don't copy untill this perticular enzyme is copied + # case: This enzyme_cplx might exist in modelA if enzyme exist + # which will automatically copie's the pool + copied = False + + if pool.parent.className not in ["Enz","ZombieEnz"]: + poolcopied = moose.copy(pool,path) + copied = True + # deleting function and enzyme which gets copied if exist under pool + # This is done to ensure daggling function / enzyme not copied. + funclist = [] + for types in ['setConc','setN','increment']: + + funclist.extend(moose.element(poolcopied).neighbors[types]) + for fl in funclist: + moose.delete(fl) + enzlist = moose.element(poolcopied).neighbors['reac'] + for el in list(set(enzlist)): + moose.delete(el.path) + return copied + +def updatePoolList(comptAdict): + for key,value in list(comptAdict.items()): + plist = moose.wildcardFind(value.path+'/##[ISA=PoolBase]') + poolListina[key] = plist + return poolListina + +def enzymeMerge(comptA,comptB,key,poolListina): + war_msg = "" + RE_Duplicated, RE_Notcopiedyet, RE_Daggling = [], [], [] + comptApath = moose.element(comptA[key]).path + comptBpath = moose.element(comptB[key]).path + objA = moose.element(comptApath).parent.name + objB = moose.element(comptBpath).parent.name + enzyListina = moose.wildcardFind(comptApath+'/##[ISA=EnzBase]') + enzyListinb = moose.wildcardFind(comptBpath+'/##[ISA=EnzBase]') + for eb in enzyListinb: + eBsubname, eBprdname = [],[] + eBsubname = subprdList(eb,"sub") + eBprdname = subprdList(eb,"prd") + allexists, allexistp = False, False + allclean = False + + poolinAlist = poolListina[findCompartment(eb).name] + for pA in poolinAlist: + if eb.parent.name == pA.name: + eapath = eb.parent.path.replace(objB,objA) + + if not moose.exists(eapath+'/'+eb.name): + #This will take care + # -- If same enzparent name but different enzyme name + # -- or different parent/enzyme name + if eBsubname and eBprdname: + allexists = checkexist(eBsubname,objB,objA) + allexistp = checkexist(eBprdname,objB,objA) + if allexists and allexistp: + enzPool = moose.element(pA.path) + eapath = eb.parent.path.replace(objB,objA) + enz = moose.element(moose.copy(eb,moose.element(eapath))) + enzPool = enz.parent + if eb.className in ["ZombieEnz","Enz"]: + moose.connect(moose.element(enz),"enz",enzPool,"reac") + if eb.className in ["ZombieMMenz","MMenz"]: + moose.connect(enzPool,"nOut",enz,"enzDest") + connectObj(enz,eBsubname,"sub",comptA,war_msg) + connectObj(enz,eBprdname,"prd",comptA,war_msg) + allclean = True + else: + # didn't find sub or prd for this Enzyme + #print ("didn't find sub or prd for this reaction" ) + RE_Notcopiedyet.append(eb) + else: + # -- it is dagging reaction + RE_Daggling.append(eb) + #print ("This reaction \""+eb.path+"\" has no substrate/product daggling reaction are not copied") + #war_msg = war_msg+"\nThis reaction \""+eb.path+"\" has no substrate/product daggling reaction are not copied" + else: + #Same Enzyme name + # -- Same substrate and product including same volume then don't copy + # -- different substrate/product or if sub/prd's volume is different then DUPLICATE the Enzyme + allclean = False + #ea = moose.element('/'+obj+'/'+enzcompartment.name+'/'+enzparent.name+'/'+eb.name) + #ea = moose.element(pA.path+'/'+eb.name) + ea = moose.element(eb.path.replace(objB,objA)) + eAsubname = subprdList(ea,"sub") + eBsubname = subprdList(eb,"sub") + hasSamenoofsublen,hasSameS,hasSamevols = same_len_name_vol(eAsubname,eBsubname) + + eAprdname = subprdList(ea,"prd") + eBprdname = subprdList(eb,"prd") + hasSamenoofprdlen,hasSameP,hasSamevolp = same_len_name_vol(eAprdname,eBprdname) + if not all((hasSamenoofsublen,hasSameS,hasSamevols,hasSamenoofprdlen,hasSameP,hasSamevolp)): + # May be different substrate or product or volume of Sub/prd may be different, + # Duplicating the enzyme + if eBsubname and eBprdname: + allexists,allexistp = False,False + allexists = checkexist(eBsubname,objB,objA) + allexistp = checkexist(eBprdname,objB,objA) + if allexists and allexistp: + eb.name = eb.name+"_duplicated" + if eb.className in ["ZombieEnz","Enz"]: + eapath = eb.parent.path.replace(objB,objA) + enz = moose.copy(eb,moose.element(eapath)) + moose.connect(enz, 'enz', eapath, 'reac' ) + + if eb.className in ["ZombieMMenz","MMenz"]: + eapath = eb.parent.path.replace(objB,objA) + enz = moose.copy(eb.name,moose.element(eapath)) + enzinfo = moose.Annotator(enz.path+'/info') + moose.connect(moose.element(enz).parent,"nOut",moose.element(enz),"enzDest") + #moose.connect(moose.element(enz),"enz",moose.element(enz).parent,"reac") + + #moose.connect( cplxItem, 'reac', enz, 'cplx' ) + connectObj(enz,eBsubname,"sub",comptA,war_msg) + connectObj(enz,eBprdname,"prd",comptA,war_msg) + RE_Duplicated.append(enz) + allclean = True + else: + allclean = False + else: + allclean = True + + if not allclean: + # didn't find sub or prd for this enzyme + # -- it may be connected Enzyme cplx + if eBsubname and eBprdname: + RE_Notcopiedyet.append(eb) + #print ("This Enzyme \""+eb.path+"\" has no substrate/product must be connect to cplx") + #war_msg = war_msg+ "\nThis Enzyme \""+rb.path+"\" has no substrate/product must be connect to cplx" + else: + RE_Daggling.append(eb) + #print ("This enzyme \""+eb.path+"\" has no substrate/product daggling reaction are not copied") + #war_msg = war_msg+"\nThis reaction \""+eb.path+"\" has no substrate/product daggling reaction are not copied" + return RE_Duplicated,RE_Notcopiedyet,RE_Daggling + +def reacMerge(comptA,comptB,key,poolListina): + RE_Duplicated, RE_Notcopiedyet, RE_Daggling = [], [], [] + war_msg = "" + comptApath = moose.element(comptA[key]).path + comptBpath = moose.element(comptB[key]).path + objA = moose.element(comptApath).parent.name + objB = moose.element(comptBpath).parent.name + reacListina = moose.wildcardFind(comptApath+'/##[ISA=ReacBase]') + reacListinb = moose.wildcardFind(comptBpath+'/##[ISA=ReacBase]') + for rb in reacListinb: + rBsubname, rBprdname = [],[] + rBsubname = subprdList(rb,"sub") + rBprdname = subprdList(rb,"prd") + allexists, allexistp = False, False + allclean = False + + if rb.name not in [ra.name for ra in reacListina]: + # reaction name not found then copy + # And assuming that pools are copied earlier EXPECT POOL CPLX + #To be assured the it takes correct compartment name incase reaction sub's + #belongs to different compt + key = findCompartment(rb).name + if rBsubname and rBprdname: + allexists = checkexist(rBsubname,objB,objA) + allexistp = checkexist(rBprdname,objB,objA) + if allexists and allexistp: + rapath = rb.parent.path.replace(objB,objA) + reac = moose.copy(rb,moose.element(rapath)) + connectObj(reac,rBsubname,"sub",comptA,war_msg) + connectObj(reac,rBprdname,"prd",comptA,war_msg) + allclean = True + else: + # didn't find sub or prd for this reaction + # -- it may be connected Enzyme cplx + RE_Notcopiedyet.append(rb) + else: + # -- it is dagging reaction + RE_Daggling.append(rb) + #print ("This reaction \""+rb.path+"\" has no substrate/product daggling reaction are not copied") + #war_msg = war_msg+"\nThis reaction \""+rb.path+"\" has no substrate/product daggling reaction are not copied" + + else: + #Same reaction name + # -- Same substrate and product including same volume then don't copy + # -- different substrate/product or if sub/prd's volume is different then DUPLICATE the reaction + allclean = False + for ra in reacListina: + if rb.name == ra.name: + rAsubname = subprdList(ra,"sub") + rBsubname = subprdList(rb,"sub") + hasSamenoofsublen,hasSameS,hasSamevols = same_len_name_vol(rAsubname,rBsubname) + + rAprdname = subprdList(ra,"prd") + rBprdname = subprdList(rb,"prd") + hasSamenoofprdlen,hasSameP,hasSamevolp = same_len_name_vol(rAprdname,rBprdname) + if not all((hasSamenoofsublen,hasSameS,hasSamevols,hasSamenoofprdlen,hasSameP,hasSamevolp)): + # May be different substrate or product or volume of Sub/prd may be different, + # Duplicating the reaction + if rBsubname and rBprdname: + allexists,allexistp = False,False + allexists = checkexist(rBsubname,objB,objA) + allexistp = checkexist(rBprdname,objB,objA) + if allexists and allexistp: + rb.name = rb.name+"_duplicated" + #reac = moose.Reac(comptA[key].path+'/'+rb.name+"_duplicated") + rapath = rb.parent.path.replace(objB,objA) + reac = moose.copy(rb,moose.element(rapath)) + connectObj(reac,rBsubname,"sub",comptA,war_msg) + connectObj(reac,rBprdname,"prd",comptA,war_msg) + RE_Duplicated.append(reac) + allclean = True + else: + allclean = False + else: + allclean = True + + if not allclean: + # didn't find sub or prd for this reaction + # -- it may be connected Enzyme cplx + if rBsubname and rBprdname: + RE_Notcopiedyet.append(rb) + #print ("This reaction \""+rb.path+"\" has no substrate/product must be connect to cplx") + #war_msg = war_msg+ "\nThis reaction \""+rb.path+"\" has no substrate/product must be connect to cplx" + else: + RE_Daggling.append(rb) + #print ("This reaction \""+rb.path+"\" has no substrate/product daggling reaction are not copied") + #war_msg = war_msg+"\nThis reaction \""+rb.path+"\" has no substrate/product daggling reaction are not copied" + + return RE_Duplicated,RE_Notcopiedyet,RE_Daggling + +def subprdList(reac,subprd): + #print "Reac ",reac + rtype = moose.element(reac).neighbors[subprd] + rname = [] + for rs in rtype: + rname.append(moose.element(rs)) + return rname + +def same_len_name_vol(rA,rB): + uaS = set(rA) + ubS = set(rB) + aS = set([uas.name for uas in uaS]) + bS = set([ubs.name for ubs in ubS]) + + hassameLen = False + hassameSP = False + hassamevol = False + hassamevollist = [] + if (len(rA) == len(rB) ): + hassameLen = True + if not (len (aS.union(bS) - aS.intersection(bS))): + hassameSP = True + if rB and rA: + rAdict = dict( [ (i.name,i) for i in (rA) ] ) + rBdict = dict( [ (i.name,i) for i in (rB) ] ) + + for key,bpath in rBdict.items(): + apath = rAdict[key] + comptA = moose.element(findCompartment(apath)) + comptB = moose.element(findCompartment(bpath)) + if not abs(comptA.volume -comptB.volume): + hassamevollist.append(True) + else: + hassamevollist.append(False) + + if len(set(hassamevollist))==1: + for x in set(hassamevollist): + hassamevol = x + + return ( hassameLen,hassameSP,hassamevol) + +def connectObj(reac,spList,spType,comptA,war_msg): + #It should not come here unless the sub/prd is connected to enzyme cplx pool + allclean = False + for rsp in spList: + for akey in list(poolListina[findCompartment(rsp).name]): + if rsp.name == akey.name: + if moose.exists(akey.path): + moose.connect(moose.element(reac), spType, moose.element(akey), 'reac', 'OneToOne') + allclean = True + else: + #It should not come here unless the sub/prd is connected to enzyme cplx pool + #print ("This pool \""+rsp.name+"\" doesnot exists in this "+comptName+" compartment to connect to this reaction \""+reac.name+"\"") + #war_msg = war_msg+ "This pool \""+rsp.name+"\" doesnot exists in this "+comptName+" compartment to connect to this reaction \""+reac.name+"\"" + allclean = False + return allclean + +def checkexist(spList,objB,objA): + allexistL = [] + allexist = False + for rsp in spList: + found = False + rspPath = rsp.path.replace(objB,objA) + if moose.exists(rspPath): + found = True + allexistL.append(found) + + if len(set(allexistL))==1: + for x in set(allexistL): + allexist = x + + return allexist + +def findCompartment(element): + while not mooseIsInstance(element,["CubeMesh","CyclMesh"]): + element = element.parent + return element + +def mooseIsInstance(element, classNames): + return moose.element(element).__class__.__name__ in classNames + + +if __name__ == "__main__": + + modelA = '/home/harsha/genesis_files/gfile/acc92.g' + modelB = '/home/harsha/genesis_files/gfile/acc50.g' + mergered = mergeChemModel(modelA,modelB) diff --git a/moose-core/python/moose/merge/mtypes.py b/moose-core/python/moose/merge/mtypes.py new file mode 100644 index 0000000000000000000000000000000000000000..9d04c5a2456ada2275dad39ed5fb2574d26abce4 --- /dev/null +++ b/moose-core/python/moose/merge/mtypes.py @@ -0,0 +1,245 @@ +# mtypes.py --- +# +# Filename: mtypes.py +# Description: +# Author: +# Maintainer: +# Created: Fri Feb 8 11:29:36 2013 (+0530) +# Version: +# Last-Updated: Tue Mar 1 02:52:35 2016 (-0500) +# By: subha +# Update #: 182 +# URL: +# Keywords: +# Compatibility: +# +# + +# Commentary: +# +# Utility to detect the model type in a file +# +# + +# Change log: +# +# +# +# +# This program is free software; you can redistribute it and/or +# modify it under the terms of the GNU General Public License as +# published by the Free Software Foundation; either version 3, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +# General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program; see the file COPYING. If not, write to +# the Free Software Foundation, Inc., 51 Franklin Street, Fifth +# Floor, Boston, MA 02110-1301, USA. +# +# + +# Code: + +from __future__ import print_function +import re +import moose + +##!!!!!!!!!!!!!!!!!!!!!!!!!!!! +## !! This was stolen from +## http://eli.thegreenplace.net/2011/10/19/perls-guess-if-file-is-text-or-binary-implemented-in-python/#id2 + +import sys +PY3 = sys.version_info[0] >= 3 + +# A function that takes an integer in the 8-bit range and returns +# a single-character byte object in py3 / a single-character string +# in py2. +# +int2byte = (lambda x: bytes((x,))) if PY3 else chr + +_text_characters = ( + b''.join(int2byte(i) for i in range(32, 127)) + + b'\n\r\t\f\b') + +def istextfile(fileobj, blocksize=512): + """ Uses heuristics to guess whether the given file is text or binary, + by reading a single block of bytes from the file. + If more than 30% of the chars in the block are non-text, or there + are NUL ('\x00') bytes in the block, assume this is a binary file. + + - Originally written by Eli Bendersky from the algorithm used + in Perl: + http://eli.thegreenplace.net/2011/10/19/perls-guess-if-file-is-text-or-binary-implemented-in-python/#id2 + """ + block = fileobj.read(blocksize) + if b'\x00' in block: + # Files with null bytes are binary + return False + elif not block: + # An empty file is considered a valid text file + return True + + # Use translate's 'deletechars' argument to efficiently remove all + # occurrences of _text_characters from the block + nontext = block.translate(None, _text_characters) + return float(len(nontext)) / len(block) <= 0.30 + +## Eli Bendersky till here !! +##!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + +def getType(filename, mode='t'): + """Returns the type of the model in file `filename`. Returns None + if type is not known. + + mode: 'b' for binary, 't' for text. Not used currently. + + """ + mtype = None + msubtype = None + if mode == 't': + for typename, typefunc in list(typeChecks.items()): + if typefunc(filename): + return typename + return None + +def getSubtype(filename, typename): + """Returns what subtype of the specified `typename` is the model file. + None if could not resolve the subtype. + """ + for subtype in subtypes[typename]: + subtypeFunc = subtypeChecks['%s/%s' % (typename, subtype)] + if subtypeFunc(filename): + return subtype + return '' + +# Dictionary of model description types and functions to detect them. +# +# Fri Feb 8 11:07:41 IST 2013 - as of now we only recognize GENESIS .g +# XML .xml extensions +# +# Fri Feb 22 16:29:43 IST 2013 - added .cspace + + + +isCSPACE = lambda x: x.lower().endswith('.cspace') +isGENESIS = lambda x: x.lower().endswith('.g') or x.endswith('.p') +isXML = lambda x: x.lower().endswith('.xml') + +isProto = lambda x: x.lower().endswith('.p') + +import xml.dom.minidom as md + +def isNeuroML(filename): + """Check if a model is in neuroml format. An xml document is + considered a neuroml if the top level element is either + 'networkml', morphml', 'channelml' or 'neuroml'. + + """ + doc = md.parse(filename) + for child in doc.childNodes: + print(child.nodeName, child.nodeType == child.ELEMENT_NODE) + if child.nodeType == child.ELEMENT_NODE and \ + (child.nodeName == 'networkml' or \ + child.nodeName == 'morphml' or \ + child.nodeName == 'channelml'or \ + child.nodeName == 'neuroml'): + return True + return False + +def isSBML(filename): + """Check model in `filename` is in SBML format.""" + doc = md.parse(filename) + for child in doc.childNodes: + if child.nodeType == child.ELEMENT_NODE and \ + child.nodeName == 'sbml': + return True + return False + +def isKKIT(filename): + """Check if `filename` is a GENESIS/KINETIKIT file. + + """ + pattern = re.compile('include\s+kkit') # KKIT files must have "include kkit" statement somewhere + with open(filename, 'r') as infile: + while True: + sentence = '' + contd = False # Flags if we are inside a multi-line entry + line = infile.readline() + if not line: # End of file + return False + line = line.strip() + # print 'read:', line + if line.find('//') == 0: # skip c++ style comment lines + # print 'c++ comment' + continue + # Skip C style multi-line comments + comment_start = line.find('/*') + if comment_start >= 0: + sentence = line[:comment_start] + # print 'cstart', comment_start, sentence + while comment_start >= 0 and line: + # print '#', line + comment_end = line.find('*/') + if comment_end >= 0: + comment_start = -1; + line = line[comment_end+2:] # keep the rest of the line + break + line = infile.readline() + if line: + line = line.strip() + # Keep extending the sentence with next line if current + # line ends with continuation token '\' + if line: + contd = line.endswith('\\') + while line and contd: + sentence += ' ' + line[:-1] + line = infile.readline() + if line: + line = line.strip() + contd = line.endswith('\\') + # if contd turned false, the last line came out of the + # while loop unprocessed + if line: + sentence += ' ' + line + # print '>', sentence + if re.search(pattern, sentence): + return True + return False + +# Mapping model types to functions to check them. These functions +# should take the filename as the single argument and return True or +# False. Define new functions for new model type checks and add them +# here. +typeChecks = { + 'cspace': isCSPACE, + 'genesis': isGENESIS, + 'xml': isXML, +} + +# These are "type/subtype": function maps for checking model subtype +# New model types with subtypes should be added here with +# corresponding checker function. +subtypeChecks = { + 'genesis/kkit': isKKIT, + 'genesis/proto': isProto, + 'xml/neuroml': isNeuroML, + 'xml/sbml': isSBML, +} + +# Mapping types to list of subtypes. +subtypes = { + 'cspace': [], + 'genesis': ['kkit', 'proto'], + 'xml': ['neuroml', 'sbml'], +} + + + +# +# mtypes.py ends here diff --git a/moose-core/python/moose/moose.py b/moose-core/python/moose/moose.py index edc9954ede60f22beac919c43c91e22210dd60b6..cda22e8f7b15c42b27cd0086991260d719c503d9 100644 --- a/moose-core/python/moose/moose.py +++ b/moose-core/python/moose/moose.py @@ -33,7 +33,7 @@ except Exception as e: print('\tError was %s' % e) genesisSupport_ = False -from . import add_Delete_ChemicalSolver +import chemUtil.add_Delete_ChemicalSolver sequence_types = ['vector<double>', 'vector<int>', @@ -106,7 +106,7 @@ def mooseWriteSBML(modelpath, filenpath, sceneitems={}): return SBML.writeSBML.mooseWriteSBML(modelpath, filepath, sceneitems) -def mooseWriteKkit(modelpath, filepath): +def mooseWriteKkit(modelpath, filepath,sceneitems={}): """Writes loded model under modelpath to a file in Kkit format. keyword arguments:\n @@ -119,7 +119,7 @@ def mooseWriteKkit(modelpath, filepath): print('GENESIS(kkit) support was not loaded') return None - return moose.genesis.writeKkit.mooseWiteKkit(modelpath, filepath) + return genesis.writeKkit.mooseWriteKkit(modelpath, filepath,sceneitems) def moosedeleteChemSolver(modelpath): @@ -128,7 +128,7 @@ def moosedeleteChemSolver(modelpath): this should be followed by mooseaddChemSolver for add solvers on to compartment to simulate else default is Exponential Euler (ee) """ - return add_Delete_ChemicalSolver.moosedeleteChemSolver(modelpath) + return chemUtil.add_Delete_ChemicalSolver.moosedeleteChemSolver(modelpath) def mooseaddChemSolver(modelpath, solver): @@ -142,7 +142,7 @@ def mooseaddChemSolver(modelpath, solver): "Runge Kutta" ("gsl") """ - return add_Delete_ChemicalSolver.mooseaddChemSolver(modelpath, solver) + return chemUtil.add_Delete_ChemicalSolver.mooseaddChemSolver(modelpath, solver) ################################################################ # Wrappers for global functions