SBMLImporter Class Reference

#include <SBMLImporter.h>

Collaboration diagram for SBMLImporter:

Public Member Functions
void	deleteCopasiModel ()
std::string	findIdInASTTree (const ASTNode *pMath, const std::set< std::string > &reactionIds)
std::string	findIdInASTTree (const ASTNode *pMath, const std::map< std::string, double > &reactionIds)
bool	getImportCOPASIMIRIAM () const
CProcessReport *	getImportHandlerAddr ()
CModel *	parseSBML (const std::string &sbmlDocumentText, CFunctionDB *funDB, SBMLDocument *&pSBMLDocument, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, CListOfLayouts *&prLol, CCopasiDataModel *pDataModel)
CModel *	readSBML (std::string filename, CFunctionDB *funDB, SBMLDocument *&pSBMLDocument, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, CListOfLayouts *&prLol, CCopasiDataModel *pDataModel)
void	restoreFunctionDB ()
	SBMLImporter ()
void	setImportCOPASIMIRIAM (bool import)
void	setImportHandler (CProcessReport *pHandler)
	~SBMLImporter ()

Static Public Member Functions
static bool	areEqualFunctions (const CFunction *pFun, const CFunction *pFun2)
static bool	areEqualSubtrees (const CEvaluationNode *pNode1, const CEvaluationNode *pNode2)
static bool	areSBMLUnitDefinitionsIdentical (const UnitDefinition *pUdef1, const UnitDefinition *pUdef2)
static CEvaluationNode *	divideByObject (const CEvaluationNode *pOrigNode, const CCopasiObject *pObject)
static UnitDefinition *	getSBMLUnitDefinitionForId (const std::string &unitId, const Model *pSBMLModel)
static bool	importNotes (CAnnotation *pAnno, const SBase *pSBase)

Protected Member Functions
void	applyStoichiometricExpressions (std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Model *pSBMLModel)
void	areRulesUnique (const Model *copasiMode)
void	checkElementUnits (const Model *pSBMLModel, CModel *pCopasiModel, int level, int version)
void	checkRuleMathConsistency (const Rule *pRule, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
bool	containsVolume (const ASTNode *pNode, const std::string &compartmentCN)
std::map< std::string, ASTNode * >	createBVarMap (const ASTNode *uDefFunction, const ASTNode *function)
CCompartment *	createCCompartmentFromCompartment (const Compartment *sbmlComp, CModel *copasiModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, const Model *pSBMLModel)
CFunction *	createCFunctionFromFunctionDefinition (const FunctionDefinition *sbmlFunction, CFunctionDB *pTmpFunctionDB, Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
CFunction *	createCFunctionFromFunctionTree (const FunctionDefinition *pSBMLFunction, Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
CMetab *	createCMetabFromSpecies (const Species *sbmlSpecies, CModel *copasiModel, CCompartment *copasiCompartment, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, const Model *pSBMLModel)
CModel *	createCModelFromSBMLDocument (SBMLDocument *doc, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
CModelValue *	createCModelValueFromParameter (const Parameter *sbmlParameter, CModel *copasiModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
CReaction *	createCReactionFromReaction (Reaction *sbmlReaction, Model *sbmlModel, CModel *cmodel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, CFunctionDB *pTmpFunctionDB)
void	createDelayFunctionDefinition ()
CEvaluationTree *	createExpressionFromFunction (const CFunction *pFun, const std::vector< std::vector< std::string > > &functionArgumentCNs)
void	createHasOnlySubstanceUnitFactor (Model *pSBMLModel, double factor, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
void	doMapping (CReaction *pCopasiReaction, const CEvaluationNodeCall *pCallNode)
void	find_local_parameters_in_delay (ASTNode *pNode, Reaction *pSBMLReaction, Model *pModel, std::map< std::string, std::string > &localReplacementMap, const std::set< std::string > &localIds, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
void	findAvogadroConstant (Model *pSBMLModel, double factor)
CFunction *	findCorrespondingFunction (const CExpression *pExpression, const CReaction *reaction)
void	findFunctionCalls (const CEvaluationNode *pNode, std::set< std::string > &functionNameSet)
const FunctionDefinition *	getFunctionDefinitionForName (const std::string name, const Model *model)
void	getIdsFromNode (const ASTNode *pNode, std::set< std::string > &idSet)
std::pair< CModel::AreaUnit, bool >	handleAreaUnit (const UnitDefinition *uDef)
std::pair< CModel::LengthUnit, bool >	handleLengthUnit (const UnitDefinition *uDef)
std::pair < CModel::QuantityUnit, bool >	handleSubstanceUnit (const UnitDefinition *uDef)
std::pair< CModel::TimeUnit, bool >	handleTimeUnit (const UnitDefinition *uDef)
std::pair< CModel::VolumeUnit, bool >	handleVolumeUnit (const UnitDefinition *uDef)
void	importEvent (const Event *pEvent, Model *pSBMLModel, CModel *pCopasiModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
void	importEvents (Model *pSBMLModel, CModel *pCopasiModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
CFunctionDB *	importFunctionDefinitions (Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
void	importInitialAssignments (Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlMap, const CModel *pCOPASIModel)
bool	importMIRIAM (const SBase *pSBMLObject, CCopasiObject *pCOPASIObject)
void	importRule (const Rule *rule, CModelEntity::Status ruleType, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Model *pSBMLModel)
void	importRuleForModelEntity (const Rule *rule, CModelEntity *pMV, CModelEntity::Status ruleType, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Model *pSBMLModel)
void	importSBMLRule (const Rule *sbmlRule, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Model *pSBMLModel)
CCopasiObject *	isConstantFlux (const CEvaluationNode *pRoot, CModel *pModel, CFunctionDB *pFunctionDB)
bool	isDelayFunctionUsed (ConverterASTNode *pNode)
std::vector < CEvaluationNodeObject * > *	isMassAction (const CEvaluationTree *pTree, const CChemEq &chemicalEquation, const CEvaluationNodeCall *pCallNode=NULL)
std::vector < CEvaluationNodeObject * > *	isMassActionExpression (const CEvaluationNode *pRootNode, const CChemEq &chemicalEquation)
std::vector < CEvaluationNodeObject * > *	isMassActionFunction (const CFunction *pFun, const CChemEq &chemicalEquation, const std::vector< std::vector< std::string > > &functionArgumentCNs)
ConverterASTNode *	isMultipliedByVolume (const ASTNode *node, const std::string &compartmentSBMLId)
bool	isSimpleFunctionCall (const CEvaluationNode *pRootNode)
bool	isStochasticModel (const Model *pSBMLModel)
void	multiplySubstanceOnlySpeciesByVolume (ConverterASTNode *pNode)
void	preprocessNode (ConverterASTNode *pNode, Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Reaction *pSBMLReaction=NULL)
bool	removeUnusedFunctions (CFunctionDB *pTmpFunctionDB, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
void	renameMassActionParameters (CEvaluationNodeCall *pCallNode)
void	replace_delay_nodes (ConverterASTNode *pNode, Model *pModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Reaction *pSBMLReaction, std::map< std::string, std::string > &localReplacementMap)
void	replace_name_nodes (ASTNode *pNode, const std::map< std::string, std::string > &replacementMap)
void	replace_time_with_initial_time (ASTNode *pNode, const CModel *pCOPASIModel)
void	replaceAmountReferences (ConverterASTNode *pNode, Model *pSBMLModel, double factor, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
void	replaceCallNodeNames (ASTNode *pNode)
void	replaceLog (ConverterASTNode *sourceNode)
void	replaceObjectNames (ASTNode *pNode, const std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, bool initialExpression=false)
void	replacePowerFunctionNodes (ASTNode *node)
void	replaceRoot (ConverterASTNode *sourceNode)
void	replaceTimeAndAvogadroNodeNames (ASTNode *pNode)
void	replaceTimeDependentFunctionCalls (ASTNode *root)
bool	replaceTimeNodesInFunctionDefinition (ASTNode *root, std::string newNodeName)
bool	sbmlId2CopasiCN (ASTNode *pNode, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, CCopasiParameterGroup &pParamGroup)
void	separateProductArguments (const CEvaluationNode *pRootNode, std::vector< const CEvaluationNode * > &arguments)
void	setCorrectUsage (CReaction *pCopasiReaction, const CEvaluationNodeCall *pCallNode)
bool	setInitialValues (CModel *pModel, const std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
CEvaluationNode *	variables2objects (const CEvaluationNode *pOrigNode, const std::map< std::string, std::string > &replacementMap)

Static Protected Member Functions
static bool	areApproximatelyEqual (const double &x, const double &y, const double &t=1e-9)
static Unit *	convertSBMLCubicmetresToLitres (const Unit *pU)
static void	findDirectDependencies (const FunctionDefinition *pFunDef, std::map< const FunctionDefinition *, std::set< std::string > > &dependencies)
static void	findDirectDependencies (const ASTNode *pNode, std::set< std::string > &dependencies)
static void	normalizeSBMLUnit (Unit *pU)
static C_FLOAT64	round (const C_FLOAT64 &x)

Protected Attributes
CFunctionDB *	functionDB
bool	mAssignmentToSpeciesReferenceFound
bool	mAvogadroCreated
bool	mAvogadroSet
bool	mDelayFound
std::map< std::string, std::string >	mDelayNodeMap
std::set< std::string >	mDivisionByCompartmentReactions
std::set< std::string >	mExplicitelyTimeDependentFunctionDefinitions
std::set< std::string >	mFastReactions
std::map< std::string, std::string >	mFunctionNameMapping
size_t	mhImportStep
std::vector< std::string >	mIgnoredParameterUnits
std::set< unsigned int >	mIgnoredSBMLMessages
bool	mImportCOPASIMIRIAM
unsigned C_INT32	mImportStep
bool	mIncompleteModel
std::map< std::string, std::string >	mKnownCustomUserDefinedFunctions
std::map< std::string, std::string >	mKnownInitalValues
unsigned int	mLevel
unsigned int	mOriginalLevel
CModel *	mpCopasiModel
CCopasiDataModel *	mpDataModel
CProcessReport *	mpImportHandler
std::set< const Parameter * >	mPotentialAvogadroNumbers
std::set< std::string >	mReactionsWithReplacedLocalParameters
std::map< const ASTNode *, CChemEqElement * >	mStoichiometricExpressionMap
std::map< Species *, Compartment * >	mSubstanceOnlySpecies
unsigned C_INT32	mTotalSteps
bool	mUnitOnNumberFound
bool	mUnsupportedAssignmentRuleFound
bool	mUnsupportedRateRuleFound
bool	mUnsupportedRuleFound
std::set< std::string >	mUsedFunctions
std::set< std::string >	mUsedSBMLIds
bool	mUsedSBMLIdsPopulated
unsigned int	mVersion
std::map< CFunction *, std::string >	sbmlIdMap
std::map< std::string, CMetab * >	speciesMap

Detailed Description

Definition at line 44 of file SBMLImporter.h.

Constructor & Destructor Documentation

SBMLImporter::SBMLImporter

(

)

Constructor that initializes speciesMap and the FunctionDB object

Definition at line 3045 of file SBMLImporter.cpp.

References functionDB, mAssignmentToSpeciesReferenceFound, mAvogadroCreated, mDelayFound, mIgnoredSBMLMessages, mIncompleteModel, mpImportHandler, mUnitOnNumberFound, mUnsupportedAssignmentRuleFound, mUnsupportedRateRuleFound, mUnsupportedRuleFound, and speciesMap.

                          :
  mIgnoredSBMLMessages(),
  speciesMap(),
  functionDB(NULL),
  mIncompleteModel(false),
  mUnsupportedRuleFound(false),
  mUnsupportedRateRuleFound(false),
  mUnsupportedAssignmentRuleFound(false),
  mUnitOnNumberFound(false),
  mAssignmentToSpeciesReferenceFound(false),
  mLevel(0),
  mOriginalLevel(0),
  mVersion(0),
  sbmlIdMap(),
  mUsedFunctions(),
  mpDataModel(NULL),
  mpCopasiModel(NULL),
  mFunctionNameMapping(),
  mDivisionByCompartmentReactions(),
  mpImportHandler(NULL),
  mImportStep(0),
  mhImportStep(C_INVALID_INDEX),
  mTotalSteps(0),
  mSubstanceOnlySpecies(),
  mFastReactions(),
  mReactionsWithReplacedLocalParameters(),
  mExplicitelyTimeDependentFunctionDefinitions(),
  mIgnoredParameterUnits(),
  mStoichiometricExpressionMap(),
  mDelayFound(false),
  mPotentialAvogadroNumbers(),
  mAvogadroCreated(false),
  mImportCOPASIMIRIAM(true),
  mDelayNodeMap(),
  mUsedSBMLIds(),
  mUsedSBMLIdsPopulated(false),
  mAvogadroSet(false),
  mKnownCustomUserDefinedFunctions(),
#if LIBSBML_VERSION >= 40100
  mpModelConversionFactor(NULL),
  mChemEqElementSpeciesIdMap(),
  mSpeciesConversionParameterMap(),
  mSBMLIdModelValueMap(),
  mSBMLSpeciesReferenceIds(),
  mRateRuleForSpeciesReferenceIgnored(false),
  mEventAssignmentForSpeciesReferenceIgnored(false),
  mConversionFactorFound(false),
# if LIBSBML_VERSION >= 40200
  mEventPrioritiesIgnored(false),
  mInitialTriggerValues(false),
  mNonPersistentTriggerFound(false)
# endif // LIBSBML_VERSION >= 40200
#endif // LIBSBML_VERSION >= 40100
{
  this->speciesMap = std::map<std::string, CMetab*>();
  this->functionDB = NULL;
  this->mIncompleteModel = false;
  this->mUnsupportedRuleFound = false;
  this->mUnsupportedRateRuleFound = false;
  this->mUnsupportedAssignmentRuleFound = false;
  this->mUnitOnNumberFound = false;
  this->mAssignmentToSpeciesReferenceFound = false;
  this->mpImportHandler = NULL;
  this->mDelayFound = false;
  this->mAvogadroCreated = false;
  // these data structures are used to handle the new conversion factores in
  // SBML L3 models and references to species references
  this->mpModelConversionFactor = NULL;
  this->mChemEqElementSpeciesIdMap.clear();
  this->mSpeciesConversionParameterMap.clear();
  this->mSBMLIdModelValueMap.clear();
  this->mRateRuleForSpeciesReferenceIgnored = false;
  this->mEventAssignmentForSpeciesReferenceIgnored = false;
  this->mConversionFactorFound = false;
  this->mEventPrioritiesIgnored = false;
  this->mInitialTriggerValues = false;
  this->mNonPersistentTriggerFound = false;

  this->mIgnoredSBMLMessages.insert(10501);
  this->mIgnoredSBMLMessages.insert(10512);
  this->mIgnoredSBMLMessages.insert(10513);
  this->mIgnoredSBMLMessages.insert(10522);
  this->mIgnoredSBMLMessages.insert(10533);
  this->mIgnoredSBMLMessages.insert(10541);
  this->mIgnoredSBMLMessages.insert(10551);
  this->mIgnoredSBMLMessages.insert(10562);
  this->mIgnoredSBMLMessages.insert(80701);
  this->mIgnoredSBMLMessages.insert(99505);
}

SBMLImporter::~SBMLImporter

(

)

Destructor that does nothing.

Definition at line 3138 of file SBMLImporter.cpp.

{}

Member Function Documentation

void SBMLImporter::applyStoichiometricExpressions ( std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap, Model *  pSBMLModel  )

protected

This method evaluates all stoichiometric expressions and sets them as constants on the CChemEqElement.

Definition at line 8980 of file SBMLImporter.cpp.

References CChemEq::addMetabolite(), CCopasiVector< T >::begin(), CExpression::calcValue(), CExpression::compile(), CCopasiVector< T >::end(), fatalError, CChemEqElement::getMetabolite(), CChemEqElement::getMetaboliteKey(), CChemEqElement::getMultiplicity(), CCopasiObject::getObjectParent(), CEvaluationTree::getRoot(), CModelEntity::getSBMLId(), CReaction::getSBMLId(), CChemEq::getSubstrates(), MCSBML, mpCopasiModel, mpDataModel, mStoichiometricExpressionMap, preprocessNode(), CChemEq::PRODUCT, replaceObjectNames(), CEvaluationTree::setTree(), CChemEq::SUBSTRATE, and CCopasiMessage::WARNING.

Referenced by createCModelFromSBMLDocument().

{
  bool warningDone = false;
  std::map<const ASTNode*, CChemEqElement* >::iterator it = this->mStoichiometricExpressionMap.begin(), end = this->mStoichiometricExpressionMap.end();
  std::vector<CCopasiContainer*> listOfContainers;
  listOfContainers.push_back(this->mpCopasiModel);

  while (it != end)
    {
      assert(it->second != NULL);
      CChemEqElement* pChemEqElement = it->second;
      ConverterASTNode* pNode = new ConverterASTNode(*it->first);
      this->preprocessNode(pNode, pSBMLModel, copasi2sbmlmap);
      this->replaceObjectNames(pNode, copasi2sbmlmap, true);
      CExpression* pExpr = new CExpression("", mpDataModel);
      pExpr->setTree(*pNode);
      pExpr->compile(listOfContainers);
      delete pNode;

      if (pExpr->getRoot() == NULL)
        {
          const CReaction* pR = dynamic_cast<const CReaction*>(pChemEqElement->getObjectParent()->getObjectParent()->getObjectParent());
          std::string id = pChemEqElement->getMetabolite()->getSBMLId();
          // create an error message that some stoichiometric expression
          // could not be evaluated that the value for the stoichiometry has
          // been set to 1.0
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 65, id.c_str(), pR->getSBMLId().c_str());
        }
      else
        {
          double value = pExpr->calcValue();
          value -= pChemEqElement->getMultiplicity();
          // find out if the metabolite is a substrate or a product
          delete pExpr;
          CChemEq* pChemEq = dynamic_cast<CChemEq*>(pChemEqElement->getObjectParent()->getObjectParent());
          assert(pChemEq != NULL);

          if (pChemEq != NULL)
            {
              CCopasiVector < CChemEqElement >::const_iterator iit = pChemEq->getSubstrates().begin(), iendit = pChemEq->getSubstrates().end();

              while (iit != iendit)
                {
                  if ((*iit) == pChemEqElement)
                    {
                      break;
                    }

                  ++iit;
                }

              if (iit != iendit)
                {
                  pChemEq->addMetabolite(pChemEqElement->getMetaboliteKey(), value, CChemEq::SUBSTRATE);
                }
              else
                {
                  pChemEq->addMetabolite(pChemEqElement->getMetaboliteKey(), value, CChemEq::PRODUCT);
                }

              // give a warning that an stoichiometric expression has been
              // converted into a constant
              if (!warningDone && !this->mStoichiometricExpressionMap.empty())
                {
                  CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 64);
                  warningDone = true;
                }
            }
          else
            {
              fatalError();
            }
        }

      ++it;
    }
}

bool SBMLImporter::areApproximatelyEqual ( const double &  x, const double &  y, const double &  t = 1e-9  )

staticprotected

Definition at line 177 of file SBMLImporter.cpp.

References min.

Referenced by areSBMLUnitDefinitionsIdentical(), findAvogadroConstant(), handleAreaUnit(), handleLengthUnit(), handleSubstanceUnit(), handleTimeUnit(), handleVolumeUnit(), isMassActionExpression(), and replaceAmountReferences().

{
  double Scale =
    (fabs(x) + fabs(y)) * t;

  // Avoid underflow
  if (Scale < 100.0 * std::numeric_limits< C_FLOAT64 >::min())
    return true;

  return 2 * fabs(x - y) < Scale;
}

bool SBMLImporter::areEqualFunctions ( const CFunction *  pFun, const CFunction *  pFun2  )

static

Definition at line 5121 of file SBMLImporter.cpp.

References areEqualSubtrees(), CCopasiObject::getObjectName(), CEvaluationTree::getRoot(), CFunction::getVariables(), and CFunctionParameters::size().

Referenced by createCFunctionFromFunctionDefinition(), findCorrespondingFunction(), and CReaction::setFunctionFromExpressionTree().

{
  bool result = true;
  const CFunctionParameters& funParams1 = pFun->getVariables();
  const CFunctionParameters& funParams2 = pFun2->getVariables();

  if (funParams1.size() == funParams2.size())
    {
      size_t i, iMax = funParams1.size();

      for (i = 0; i < iMax; ++i)
        {
          const CFunctionParameter* pFunParam1 = funParams1[i];
          const CFunctionParameter* pFunParam2 = funParams2[i];

          if (pFunParam1->getObjectName() != pFunParam2->getObjectName())
            {
              result = false;
              break;
            }
        }

      if (result == true)
        {
          const CEvaluationNode* pNodeFun1 = static_cast<const CEvaluationNode*>(pFun->getRoot());
          const CEvaluationNode* pNodeFun2 = static_cast<const CEvaluationNode*>(pFun2->getRoot());
          result = areEqualSubtrees(pNodeFun1, pNodeFun2);
        }
    }
  else
    {
      result = false;
    }

  return result;
}

bool SBMLImporter::areEqualSubtrees ( const CEvaluationNode *  pNode1, const CEvaluationNode *  pNode2  )

static

Compares to CEvaluationNode based subtrees recursively.

Definition at line 5159 of file SBMLImporter.cpp.

References CCopasiNode< _Data >::getChild(), CCopasiNode< _Data >::getData(), CCopasiNode< _Data >::getSibling(), and CEvaluationNode::getType().

Referenced by areEqualFunctions().

{
  // TODO CRITICAL We need to use a node iterator
  bool result = ((pNode1->getType() == pNode2->getType()) && (pNode1->getData() == pNode2->getData()));
  const CEvaluationNode* pChild1 = static_cast<const CEvaluationNode*>(pNode1->getChild());
  const CEvaluationNode* pChild2 = static_cast<const CEvaluationNode*>(pNode2->getChild());

  while (result && pChild1 && pChild2)
    {
      result = areEqualSubtrees(pChild1, pChild2);
      pChild1 = static_cast<const CEvaluationNode*>(pChild1->getSibling());
      pChild2 = static_cast<const CEvaluationNode*>(pChild2->getSibling());
    }

  result = (result && !pChild1 && !pChild2);
  return result;
}

void SBMLImporter::areRulesUnique ( const Model *  copasiMode )

protected

Checks if no id is used in more than one Assignment and RateRule.

Definition at line 6662 of file SBMLImporter.cpp.

References CCopasiMessage::EXCEPTION, and MCSBML.

Referenced by createCModelFromSBMLDocument().

{
  std::set<std::string> idSet;

  // go through the rules and check that no id is used in more than one rule
  unsigned int i, iMax = sbmlModel->getNumRules();

  for (i = 0; i < iMax; ++i)
    {
      const Rule* pRule = sbmlModel->getRule(i);
      std::string id;

      switch (pRule->getTypeCode())
        {
          case SBML_ASSIGNMENT_RULE:
            id = dynamic_cast<const AssignmentRule*>(pRule)->getVariable();
            break;

          case SBML_RATE_RULE:
            id = dynamic_cast<const RateRule*>(pRule)->getVariable();
            break;

          default:
            break;
        }

      if (!id.empty())
        {
          if (!idSet.insert(id).second)
            {
              CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 35, id.c_str());
              break;
            }
        }
    }
}

bool SBMLImporter::areSBMLUnitDefinitionsIdentical ( const UnitDefinition *  pUdef1, const UnitDefinition *  pUdef2  )

static

Enhanced method to identify identical SBML unit definitions. This method uses the areIdentical method from libSBML, but if the method return false, it does some extra checks. Right now it check for example if two volumes, one given in litre and one given in cubic meters are identical.

Enhanced method to identify identical sbml unit definitions. The method first converts the unit definitions to SI units and simplifies them, only then they are compared.

Definition at line 8297 of file SBMLImporter.cpp.

References areApproximatelyEqual().

Referenced by CSBMLExporter::check_for_spatial_size_units(), checkElementUnits(), CSBMLExporter::createAreaUnit(), createCModelFromSBMLDocument(), CSBMLExporter::createLengthUnit(), CSBMLExporter::createSubstanceUnit(), CSBMLExporter::createTimeUnit(), and CSBMLExporter::createVolumeUnit().

{
  UnitDefinition* pTmpUdef1 = UnitDefinition::convertToSI(pUdef1);
  UnitDefinition::simplify(pTmpUdef1);
  UnitDefinition* pTmpUdef2 = UnitDefinition::convertToSI(pUdef2);
  UnitDefinition::simplify(pTmpUdef2);
  bool result = UnitDefinition::areIdentical(pUdef1, pUdef2);

  if (result == false)
    {
      // check if maybe everything is the same, only the multipliers are
      // somewhat off due to rounding errors
      bool newResult = true;

      if (pTmpUdef1->getNumUnits() == pTmpUdef2->getNumUnits())
        {
          UnitDefinition::reorder(pTmpUdef1);
          UnitDefinition::reorder(pTmpUdef2);
          unsigned int i = 0, iMax = pTmpUdef1->getNumUnits();

          // for COPASI mole and avogadro are the same units, so in order to compare them,
          // we replace the avogadro unit by mole
          for (i = 0; i < iMax; ++i)
            {
              if (pTmpUdef1->getUnit(i)->getKind() == UNIT_KIND_AVOGADRO)
                {
                  pTmpUdef1->getUnit(i)->setKind(UNIT_KIND_MOLE);
                }

              if (pTmpUdef2->getUnit(i)->getKind() == UNIT_KIND_AVOGADRO)
                {
                  pTmpUdef2->getUnit(i)->setKind(UNIT_KIND_MOLE);
                }
            }

          // we have to reset i to 0
          i = 0;
          const Unit *pU1, *pU2;

          while (newResult == true && i != iMax)
            {
              pU1 = pTmpUdef1->getUnit(i);
              pU2 = pTmpUdef2->getUnit(i);

              if (pU1->getKind() != pU2->getKind() ||
                  pU1->getExponent() != pU2->getExponent() ||
                  pU1->getScale() != pU2->getScale() ||
                  !areApproximatelyEqual(pU2->getMultiplier(), pU1->getMultiplier()))
                {
                  newResult = false;
                }

              ++i;
            }

          result = newResult;
        }
    }

  delete pTmpUdef1;
  delete pTmpUdef2;
  return result;
}

void SBMLImporter::checkElementUnits ( const Model *  pSBMLModel, CModel *  pCopasiModel, int  level, int  version  )

protected

Definition at line 7353 of file SBMLImporter.cpp.

References areSBMLUnitDefinitionsIdentical(), CModel::dimensionlessArea, CModel::dimensionlessLength, CModel::dimensionlessVolume, CCopasiMessage::EXCEPTION, fatalError, getSBMLUnitDefinitionForId(), handleAreaUnit(), handleLengthUnit(), handleSubstanceUnit(), handleTimeUnit(), handleVolumeUnit(), MCSBML, CModel::setAreaUnit(), CModel::setLengthUnit(), CModel::setTimeUnit(), CModel::setVolumeUnit(), and CCopasiMessage::WARNING.

Referenced by createCModelFromSBMLDocument().

{
  unsigned int i, iMax = pSBMLModel->getNumCompartments();
  std::vector<std::string> nonDefaultCompartments;
  std::vector<std::string> nonDefaultSpecies;
  std::vector<std::string> nonDefaultKineticTime;
  std::vector<std::string> nonDefaultKineticSubstance;
  std::vector<std::string> nonDefaultEventTime;
  const Compartment* pCompartment;
  const Species* pSpecies;
  const Reaction* pReaction;
  const KineticLaw* pKineticLaw;
  UnitDefinition* pDimensionlessUnits = getSBMLUnitDefinitionForId("dimensionless", pSBMLModel);
  UnitDefinition* pLengthUnits = getSBMLUnitDefinitionForId("length", pSBMLModel);
  UnitDefinition* pAreaUnits = getSBMLUnitDefinitionForId("area", pSBMLModel);
  UnitDefinition* pVolumeUnits = getSBMLUnitDefinitionForId("volume", pSBMLModel);
  UnitDefinition* pTimeUnits = getSBMLUnitDefinitionForId("time", pSBMLModel);
  UnitDefinition* pSubstanceUnits = getSBMLUnitDefinitionForId("substance", pSBMLModel);
  std::string lastVolumeUnit = "";
  std::string lastAreaUnit = "";
  std::string lastLengthUnit = "";
  std::string lastDimensionlessUnit = "";
  bool inconsistentVolumeUnits = false;
  bool inconsistentAreaUnits = false;
  bool inconsistentLengthUnits = false;
  bool inconsistentDimensionlessUnits = false;
  bool defaultVolumeUsed = false;
  bool defaultAreaUsed = false;
  bool defaultLengthUsed = false;

  for (i = 0; i < iMax; ++i)
    {
      pCompartment = pSBMLModel->getCompartment(i);

      if (pCompartment->getSpatialDimensions() == 3)
        {
          if (pCompartment->isSetUnits())
            {
              std::string unitId = pCompartment->getUnits();
              UnitDefinition* pUdef1 = getSBMLUnitDefinitionForId(unitId, pSBMLModel);

              if (pUdef1 == NULL)
                {
                  // error message
                  CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 55, unitId.c_str(), "units", "compartment", pCompartment->getId().c_str());
                }
              else
                {
                  std::pair<CModel::VolumeUnit, bool> result = this->handleVolumeUnit(pUdef1);

                  if (result.second == false)
                    {
                      // we did not recognize the unit
                      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 79, "volume", "the 3 dimensional  compartment", pCompartment->getId().c_str());
                      delete pUdef1;
                      continue;
                    }
                }

              if (unitId != "volume" && !areSBMLUnitDefinitionsIdentical(pVolumeUnits, pUdef1))
                {
                  nonDefaultCompartments.push_back(pCompartment->getId());
                }
              else if (unitId == "volume")
                {
                  defaultVolumeUsed = true;
                }

              if (lastVolumeUnit == "")
                {
                  lastVolumeUnit = unitId;
                }
              else if (unitId != lastVolumeUnit)
                {
                  // check if the two units have identical definitions
                  UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastVolumeUnit, pSBMLModel);
                  assert(pUdef2 != NULL);

                  if (!areSBMLUnitDefinitionsIdentical(pUdef1, pUdef2))
                    {
                      inconsistentVolumeUnits = true;
                    }

                  delete pUdef2;
                }

              delete pUdef1;
            }
          else if (lastVolumeUnit == "")
            {
              lastVolumeUnit = "volume";
              defaultVolumeUsed = true;
            }
          else
            {
              defaultVolumeUsed = true;
              // compare the default volume unit to the lastVolumeUnit
              UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastVolumeUnit, pSBMLModel);
              assert(pUdef2 != NULL);

              if (!areSBMLUnitDefinitionsIdentical(pVolumeUnits, pUdef2))
                {
                  inconsistentVolumeUnits = true;
                }

              delete pUdef2;
              lastVolumeUnit = "volume";
            }
        }
      else if (pCompartment->getSpatialDimensions() == 2)
        {
          if (pCompartment->isSetUnits())
            {
              std::string unitId = pCompartment->getUnits();
              UnitDefinition* pUdef1 = getSBMLUnitDefinitionForId(unitId, pSBMLModel);

              if (pUdef1 == NULL)
                {
                  // error message
                  CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 55, unitId.c_str(), "units", "compartment", pCompartment->getId().c_str());
                }
              else
                {
                  std::pair<CModel::AreaUnit, bool> result = this->handleAreaUnit(pUdef1);

                  if (result.second == false)
                    {
                      // we did not recognize the unit
                      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 79, "area", "the 2 dimensional compartment", pCompartment->getId().c_str());
                      continue;
                    }
                }

              if (unitId != "area" && !areSBMLUnitDefinitionsIdentical(pAreaUnits, pUdef1))
                {
                  nonDefaultCompartments.push_back(pCompartment->getId());
                }
              else if (unitId == "area")
                {
                  defaultAreaUsed = true;
                }

              if (lastAreaUnit == "")
                {
                  lastAreaUnit = unitId;
                }
              else if (unitId != lastAreaUnit)
                {
                  // check if the two units have identical definitions
                  UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastAreaUnit, pSBMLModel);
                  assert(pUdef2 != NULL);

                  if (!areSBMLUnitDefinitionsIdentical(pUdef1, pUdef2))
                    {
                      inconsistentAreaUnits = true;
                    }

                  delete pUdef2;
                }

              delete pUdef1;
            }
          else if (lastAreaUnit == "")
            {
              lastAreaUnit = "area";
              defaultAreaUsed = true;
            }
          else
            {
              defaultAreaUsed = true;
              // compare the default area unit to the lastAreaUnit
              UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastAreaUnit, pSBMLModel);
              assert(pUdef2 != NULL);

              if (!areSBMLUnitDefinitionsIdentical(pAreaUnits, pUdef2))
                {
                  inconsistentAreaUnits = true;
                }

              delete pUdef2;
              lastAreaUnit = "area";
            }
        }
      else if (pCompartment->getSpatialDimensions() == 1)
        {
          if (pCompartment->isSetUnits())
            {
              std::string unitId = pCompartment->getUnits();
              UnitDefinition* pUdef1 = getSBMLUnitDefinitionForId(unitId, pSBMLModel);

              if (pUdef1 == NULL)
                {
                  // error message
                  CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 55, unitId.c_str(), "units", "compartment", pCompartment->getId().c_str());
                }
              else
                {
                  std::pair<CModel::LengthUnit, bool> result = this->handleLengthUnit(pUdef1);

                  if (result.second == false)
                    {
                      // we did not recognize the unit
                      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 79, "length", "the 1 dimensional compartment", pCompartment->getId().c_str());
                      continue;
                    }
                }

              if (unitId != "length" && !areSBMLUnitDefinitionsIdentical(pLengthUnits, pUdef1))
                {
                  nonDefaultCompartments.push_back(pCompartment->getId());
                }
              else if (unitId == "length")
                {
                  defaultLengthUsed = true;
                }

              if (lastLengthUnit == "")
                {
                  lastLengthUnit = unitId;
                }
              else if (unitId != lastLengthUnit)
                {
                  // check if the two units have identical definitions
                  UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastLengthUnit, pSBMLModel);
                  assert(pUdef2 != NULL);

                  if (!areSBMLUnitDefinitionsIdentical(pUdef1, pUdef2))
                    {
                      inconsistentLengthUnits = true;
                    }

                  delete pUdef2;
                }

              delete pUdef1;
            }
          else if (lastLengthUnit == "")
            {
              lastLengthUnit = "length";
              defaultLengthUsed = true;
            }
          else
            {
              defaultLengthUsed = true;
              // compare the default length unit to the lastLengthUnit
              UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastLengthUnit, pSBMLModel);
              assert(pUdef2 != NULL);

              if (!areSBMLUnitDefinitionsIdentical(pLengthUnits, pUdef2))
                {
                  inconsistentLengthUnits = true;
                }

              delete pUdef2;
              lastLengthUnit = "length";
            }
        }
      else if (pCompartment->getSpatialDimensions() == 0)
        {
          if (pCompartment->isSetUnits())
            {
              std::string unitId = pCompartment->getUnits();
              UnitDefinition* pUdef1 = getSBMLUnitDefinitionForId(unitId, pSBMLModel);

              if (pUdef1 == NULL)
                {
                  // error message
                  CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 55, unitId.c_str(), "units", "compartment", pCompartment->getId().c_str());
                }

              if (unitId != "dimensionless" && !areSBMLUnitDefinitionsIdentical(pDimensionlessUnits, pUdef1))
                {
                  CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 79, "dimensionless", "the 0 dimensional compartment", pCompartment->getId().c_str());
                  continue;
                }

              if (lastDimensionlessUnit == "")
                {
                  lastDimensionlessUnit = unitId;
                }
              else if (unitId != lastDimensionlessUnit)
                {
                  // check if the two units have identical definitions
                  UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastDimensionlessUnit, pSBMLModel);
                  assert(pUdef2 != NULL);

                  if (!areSBMLUnitDefinitionsIdentical(pUdef1, pUdef2))
                    {
                      inconsistentDimensionlessUnits = true;
                    }

                  delete pUdef2;
                }

              delete pUdef1;
            }
          else if (lastDimensionlessUnit == "")
            {
              lastDimensionlessUnit = "dimensionless";
            }
        }
    }

  if (!inconsistentVolumeUnits && lastVolumeUnit != "" && lastVolumeUnit != "volume")
    {
      // try to set the default volume unit to the unit defined by lastUnit
      const UnitDefinition* pUdef = SBMLImporter::getSBMLUnitDefinitionForId(lastVolumeUnit, pSBMLModel);
      assert(pUdef != NULL);
      std::pair<CModel::VolumeUnit, bool> volume = this->handleVolumeUnit(pUdef);

      if (volume.second == true)
        {
          // set the default volume unit
          pCopasiModel->setVolumeUnit(volume.first);
          delete pVolumeUnits;
          pVolumeUnits = dynamic_cast<UnitDefinition*>(pUdef->clone());
        }
      else
        {
          inconsistentVolumeUnits = true;
        }

      delete pUdef;
    }

  if (!inconsistentAreaUnits && lastAreaUnit != "" && lastAreaUnit != "area")
    {
      // try to set the default area unit to the unit defined by lastAreaUnit
      const UnitDefinition* pUdef = SBMLImporter::getSBMLUnitDefinitionForId(lastAreaUnit, pSBMLModel);
      assert(pUdef != NULL);
      std::pair<CModel::AreaUnit, bool> area = this->handleAreaUnit(pUdef);

      if (area.second == true)
        {
          // set the default area unit
          pCopasiModel->setAreaUnit(area.first);
          delete pAreaUnits;
          pAreaUnits = dynamic_cast<UnitDefinition*>(pUdef->clone());
        }
      else
        {
          inconsistentAreaUnits = true;
        }

      delete pUdef;
    }

  if (!inconsistentLengthUnits && lastLengthUnit != "" && lastLengthUnit != "length")
    {
      // try to set the default length unit to the unit defined by lastLengthUnit
      const UnitDefinition* pUdef = SBMLImporter::getSBMLUnitDefinitionForId(lastLengthUnit, pSBMLModel);
      assert(pUdef != NULL);
      std::pair<CModel::LengthUnit, bool> length = this->handleLengthUnit(pUdef);

      if (length.second == true)
        {
          // set the default length unit
          pCopasiModel->setLengthUnit(length.first);
          delete pLengthUnits;
          pLengthUnits = dynamic_cast<UnitDefinition*>(pUdef->clone());
        }
      else
        {
          inconsistentLengthUnits = true;
        }

      delete pUdef;
    }

  if (inconsistentVolumeUnits || inconsistentAreaUnits || inconsistentLengthUnits || inconsistentDimensionlessUnits)
    {
      // warn about inconsistent units and that they have been ignored and
      // report the actual units used
      // one warning for every entry in nonDefaultCompartment
      std::vector<std::string>::iterator errorIt = nonDefaultCompartments.begin(), errorEndit = nonDefaultCompartments.end();
      std::ostringstream os;

      while (errorIt != errorEndit)
        {
          os << *errorIt << ", ";
          ++errorIt;
        }

      std::string s = os.str();
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 24 , s.substr(0, s.size() - 2).c_str());
      // check if the default units have been used for any of the compartments
      // if so, the models has to use the defaults, otherwise we can just
      // choose one
      std::pair<CModel::LengthUnit, bool> length = std::pair<CModel::LengthUnit, bool>(CModel::dimensionlessLength, false);
      const UnitDefinition* pUdef = NULL;

      if (defaultLengthUsed)
        {
          pUdef = SBMLImporter::getSBMLUnitDefinitionForId("length", pSBMLModel);
          assert(pUdef != NULL);
          length = this->handleLengthUnit(pUdef);
        }
      else
        {
          if (lastLengthUnit != "")
            {
              pUdef = SBMLImporter::getSBMLUnitDefinitionForId(lastLengthUnit, pSBMLModel);
              assert(pUdef != NULL);
              length = this->handleLengthUnit(pUdef);
            }
        }

      if (length.second == true)
        {
          // set the default length unit
          pCopasiModel->setLengthUnit(length.first);
          delete pLengthUnits;
          pLengthUnits = dynamic_cast<UnitDefinition*>(pUdef->clone());
        }

      if (pUdef)
        {
          delete pUdef;
          pUdef = NULL;
        }

      std::pair<CModel::AreaUnit, bool> area = std::pair<CModel::AreaUnit, bool>(CModel::dimensionlessArea, false);

      if (defaultAreaUsed)
        {
          pUdef = SBMLImporter::getSBMLUnitDefinitionForId("area", pSBMLModel);
          assert(pUdef != NULL);
          area = this->handleAreaUnit(pUdef);
        }
      else
        {
          if (lastAreaUnit != "")
            {
              pUdef = SBMLImporter::getSBMLUnitDefinitionForId(lastAreaUnit, pSBMLModel);
              assert(pUdef != NULL);
              area = this->handleAreaUnit(pUdef);
            }
        }

      if (area.second == true)
        {
          // set the default area unit
          pCopasiModel->setAreaUnit(area.first);
          delete pAreaUnits;
          pAreaUnits = dynamic_cast<UnitDefinition*>(pUdef->clone());
        }

      if (pUdef)
        {
          delete pUdef;
          pUdef = NULL;
        }

      std::pair<CModel::VolumeUnit, bool> volume = std::pair<CModel::VolumeUnit, bool>(CModel::dimensionlessVolume, false);

      if (defaultVolumeUsed)
        {
          pUdef = SBMLImporter::getSBMLUnitDefinitionForId("volume", pSBMLModel);
          assert(pUdef != NULL);
          volume = this->handleVolumeUnit(pUdef);
        }
      else
        {
          if (lastVolumeUnit != "")
            {
              pUdef = SBMLImporter::getSBMLUnitDefinitionForId(lastVolumeUnit, pSBMLModel);
              assert(pUdef != NULL);
              volume = this->handleVolumeUnit(pUdef);
            }
        }

      if (volume.second == true)
        {
          // set the default length unit
          pCopasiModel->setVolumeUnit(volume.first);
          delete pVolumeUnits;
          pVolumeUnits = dynamic_cast<UnitDefinition*>(pUdef->clone());
        }

      if (pUdef)
        {
          delete pUdef;
          pUdef = NULL;
        }
    }

  bool inconsistentUnits = false;
  std::string lastUnit = "";

  iMax = pSBMLModel->getNumSpecies();

  for (i = 0; i < iMax; ++i)
    {
      pSpecies = pSBMLModel->getSpecies(i);

      if (level < 2 || (level == 2 && version < 3))
        {
          // check the isSetSpatialSizeUnits flag for models prior to L2V3.
          if (pSpecies->isSetSpatialSizeUnits() == true)
            {
              // check if the spatialSizeUnits is consistent with the
              // pVolumeUnits, pAreaUnits, pLengthUnits or pDimensionlessUnits
              // first we need to find the compartment
              pCompartment = pSBMLModel->getCompartment(pSpecies->getCompartment());
              assert(pCompartment != NULL);
              std::string spatialSizeUnits = pSpecies->getSpatialSizeUnits();
              UnitDefinition* pTmpUdef2 = getSBMLUnitDefinitionForId(spatialSizeUnits, pSBMLModel);

              if (pTmpUdef2 == NULL)
                {
                  // error message
                  CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 55, spatialSizeUnits.c_str(), "spatialSizeUnits", "species", pSpecies->getId().c_str());
                }

              switch (pCompartment->getSpatialDimensions())
                {
                  case 0:

                    if (!areSBMLUnitDefinitionsIdentical(pDimensionlessUnits, pTmpUdef2))
                      {
                        CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 19, pSpecies->getId().c_str());
                      }

                    break;

                  case 1:

                    if (!areSBMLUnitDefinitionsIdentical(pLengthUnits, pTmpUdef2))
                      {
                        CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 19, pSpecies->getId().c_str());
                      }

                    break;

                  case 2:

                    if (!areSBMLUnitDefinitionsIdentical(pAreaUnits, pTmpUdef2))
                      {
                        CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 19, pSpecies->getId().c_str());
                      }

                    break;

                  case 3:

                    if (!areSBMLUnitDefinitionsIdentical(pVolumeUnits, pTmpUdef2))
                      {
                        CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 19, pSpecies->getId().c_str());
                      }

                    break;

                  default:
                    fatalError();
                    break;
                }

              delete pTmpUdef2;
            }
        }

      if (pSpecies->isSetUnits())
        {
          std::string unitId = pSpecies->getUnits();
          UnitDefinition* pUdef1 = getSBMLUnitDefinitionForId(unitId, pSBMLModel);

          if (pUdef1 == NULL)
            {
              // error message
              CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 55, unitId.c_str(), "substanceUnits", "species", pSpecies->getId().c_str());
            }

          else
            {
              std::pair<CModel::QuantityUnit, bool> result = this->handleSubstanceUnit(pUdef1);

              if (result.second == false)
                {
                  // we did not recognize the unit
                  CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 79, "substance", "species", pSpecies->getId().c_str());
                  continue;
                }
            }

          if (unitId != "substance" && !areSBMLUnitDefinitionsIdentical(pSubstanceUnits, pUdef1))
            {
              nonDefaultSpecies.push_back(pSpecies->getId());
            }

          if (lastUnit == "")
            {
              lastUnit = unitId;
            }
          else if (unitId != lastUnit)
            {
              // check if the two units have identical definitions
              UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastUnit, pSBMLModel);
              assert(pUdef2 != NULL);

              if (!areSBMLUnitDefinitionsIdentical(pUdef1, pUdef2))
                {
                  inconsistentUnits = true;
                }

              delete pUdef2;
            }

          delete pUdef1;
        }
      else if (lastUnit == "")
        {
          lastUnit = "substance";
        }
    }

  bool inconsistentTimeUnits = false;
  std::string lastTimeUnits = "";
  iMax = pSBMLModel->getNumReactions();

  for (i = 0; i < iMax; ++i)
    {
      pReaction = pSBMLModel->getReaction(i);
      pKineticLaw = pReaction->getKineticLaw();

      if (pKineticLaw != NULL)
        {
          std::string unitId;

          if (pKineticLaw->isSetSubstanceUnits())
            {
              unitId = pKineticLaw->getSubstanceUnits();
              UnitDefinition* pUdef1 = getSBMLUnitDefinitionForId(unitId, pSBMLModel);

              if (pUdef1 == NULL)
                {
                  // error message
                  CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 55, unitId.c_str(), "substanceUnits", "kinetic law of the reaction", pReaction->getId().c_str());
                }

              else
                {
                  std::pair<CModel::QuantityUnit, bool> result = this->handleSubstanceUnit(pUdef1);

                  if (result.second == false)
                    {
                      // we did not recognize the unit
                      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 79, "substance", "the kinetic law in reaction", pReaction->getId().c_str());

                      continue;
                    }
                }

              if (unitId != "substance" && !areSBMLUnitDefinitionsIdentical(pSubstanceUnits, pUdef1))
                {
                  nonDefaultKineticSubstance.push_back(pReaction->getId());
                }

              if (lastUnit == "")
                {
                  lastUnit = unitId;
                }
              else if (unitId != lastUnit)
                {
                  // check if the two units have identical definitions
                  UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastUnit, pSBMLModel);
                  assert(pUdef2 != NULL);

                  if (!areSBMLUnitDefinitionsIdentical(pUdef1, pUdef2))
                    {
                      inconsistentUnits = true;
                    }

                  delete pUdef2;
                }

              delete pUdef1;
            }
          else if (lastUnit == "")
            {
              lastUnit = "substance";
            }

          if (pKineticLaw->isSetTimeUnits())
            {
              unitId = pKineticLaw->getTimeUnits();
              UnitDefinition* pUdef1 = getSBMLUnitDefinitionForId(unitId, pSBMLModel);

              if (pUdef1 == NULL)
                {
                  // error message
                  CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 55, unitId.c_str(), "timeUnits", "kinetic law of the reaction", pReaction->getId().c_str());
                }
              else
                {
                  std::pair<CModel::TimeUnit, bool> result = this->handleTimeUnit(pUdef1);

                  if (result.second == false)
                    {
                      // we did not recognize the unit
                      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 79, "time", "the kinetic law in reaction", pReaction->getId().c_str());

                      continue;
                    }
                }

              if (unitId != "time" && !areSBMLUnitDefinitionsIdentical(pTimeUnits, pUdef1))
                {
                  nonDefaultKineticTime.push_back(pReaction->getId());
                }

              if (lastTimeUnits == "")
                {
                  lastTimeUnits = unitId;
                }
              else if (unitId != lastTimeUnits)
                {
                  // check if the two units have identical definitions
                  UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastTimeUnits, pSBMLModel);
                  assert(pUdef2 != NULL);

                  if (!areSBMLUnitDefinitionsIdentical(pUdef1, pUdef2))
                    {
                      inconsistentTimeUnits = true;
                    }

                  delete pUdef2;
                }

              delete pUdef1;
            }
          else if (lastTimeUnits == "") // set the last time unit to time
            {
              lastTimeUnits = "time";
            }
        }
    }

  if (!inconsistentUnits && lastUnit != "" && lastUnit != "substance")
    {
      // we have to check if lastUnit is different from the global substance
      // unit
      // if it differs, we have to issue a warning because we can't just set
      // another global substance unit since this would change the units for
      // the kinetic laws
      UnitDefinition* pUdef = getSBMLUnitDefinitionForId(lastUnit, pSBMLModel);
      const UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId("substance", pSBMLModel);
      assert(pUdef2 != NULL);
      assert(pUdef != NULL);

      if (!areSBMLUnitDefinitionsIdentical(pUdef, pUdef2))
        {
          // warning
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 78 , "substance", "substance", "substance");
        }

      delete pUdef2;
      delete pUdef;
    }

  if (inconsistentUnits)
    {
      // warn about inconsistent units and that they have been ignored
      // one warning SBML + 25 for each species in nonDefaultSpecies
      // and one for each KineticLaw in nonDefaultKineticSubstance
      std::vector<std::string>::iterator errorIt = nonDefaultSpecies.begin(), errorEndit = nonDefaultSpecies.end();
      std::ostringstream os;

      while (errorIt != errorEndit)
        {
          os << *errorIt << ", ";
          ++errorIt;
        }

      std::string s = os.str();
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 25 , s.substr(0, s.size() - 2).c_str());
      os.str("");
      errorIt = nonDefaultKineticSubstance.begin(), errorEndit = nonDefaultKineticSubstance.end();

      while (errorIt != errorEndit)
        {
          os << *errorIt << ", ";
          ++errorIt;
        }

      s = os.str();
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 44 , s.substr(0, s.size() - 2).c_str());
    }

  if (!inconsistentTimeUnits && lastTimeUnits != "" && lastTimeUnits != "time")
    {
      // we have to check if lastUnit is different from the global substance
      // unit
      // if it differs, we have to issue a warning because we can't just set
      // another global substance unit since this would change the units for
      // the kinetic laws
      UnitDefinition* pUdef = getSBMLUnitDefinitionForId(lastTimeUnits, pSBMLModel);
      const UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId("time", pSBMLModel);
      assert(pUdef2 != NULL);
      assert(pUdef != NULL);

      if (!areSBMLUnitDefinitionsIdentical(pUdef, pUdef2))
        {
          // warning
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 78 , "time", "time", "time");
        }

      delete pUdef2;
    }

  if (inconsistentTimeUnits)
    {
      // warn about inconsistent time unit
      // one error for each entry in nonDefaultKineticTime
      std::ostringstream os;
      std::vector<std::string>::iterator errorIt = nonDefaultKineticTime.begin(), errorEndit = nonDefaultKineticTime.end();

      while (errorIt != errorEndit)
        {
          os << *errorIt << ", ";
          ++errorIt;
        }

      std::string s = os.str();
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 53, s.substr(0, s.size() - 2).c_str());
    }

  // delete the units we created
  delete pSubstanceUnits;
  delete pVolumeUnits;
  const Event* pEvent = NULL;
  inconsistentTimeUnits = false;
  iMax = pSBMLModel->getNumEvents();

  for (i = 0; i < iMax; ++i)
    {
      pEvent = pSBMLModel->getEvent(i);
      std::string unitId;

      if (pEvent->isSetTimeUnits())
        {
          unitId = pEvent->getTimeUnits();
          UnitDefinition* pUdef1 = getSBMLUnitDefinitionForId(unitId, pSBMLModel);

          if (pUdef1 == NULL)
            {
              // error message
              CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 55, unitId.c_str(), "time units", "the event", pEvent->getId().c_str());
            }

          else
            {
              std::pair<CModel::TimeUnit, bool> result = this->handleTimeUnit(pUdef1);

              if (result.second == false)
                {
                  // we did not recognize the unit
                  CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 79, "time", "the event", pEvent->getId().c_str());

                  continue;
                }
            }

          if (unitId != "time" && !areSBMLUnitDefinitionsIdentical(pTimeUnits, pUdef1))
            {
              if (pEvent->isSetId())
                {
                  nonDefaultEventTime.push_back(pEvent->getId());
                }
            }

          if (lastTimeUnits == "")
            {
              lastTimeUnits = unitId;
            }
          else if (unitId != lastTimeUnits)
            {
              // check if the two units have identical definitions
              UnitDefinition* pUdef2 = getSBMLUnitDefinitionForId(lastTimeUnits, pSBMLModel);
              assert(pUdef2 != NULL);

              if (!areSBMLUnitDefinitionsIdentical(pUdef1, pUdef2))
                {
                  inconsistentTimeUnits = true;
                }

              delete pUdef2;
            }

          delete pUdef1;
        }
      else if (lastTimeUnits == "") // set the last time unit to time
        {
          lastTimeUnits = "time";
        }
    }

  if (!inconsistentTimeUnits && lastTimeUnits != "" && lastTimeUnits != "time")
    {
      // try to set the default time units
      UnitDefinition* pUdef = getSBMLUnitDefinitionForId(lastTimeUnits, pSBMLModel);
      assert(pUdef != NULL);
      std::pair<CModel::TimeUnit, bool> time = this->handleTimeUnit(pUdef);
      delete pUdef;

      if (time.second == true)
        {
          // set the default volume unit
          pCopasiModel->setTimeUnit(time.first);
        }
      else
        {
          inconsistentTimeUnits = true;
        }
    }

  if (inconsistentTimeUnits)
    {
      // warn about inconsistent time unit
      // one error for each entry in nonDefaultKineticTime
      std::ostringstream os;
      std::vector<std::string>::iterator errorIt = nonDefaultEventTime.begin(), errorEndit = nonDefaultEventTime.end();

      while (errorIt != errorEndit)
        {
          os << *errorIt << ", ";
          ++errorIt;
        }

      std::string s = os.str();
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 53, s.substr(0, s.size() - 2).c_str());
    }

  // delete the units we created
  delete pTimeUnits;
  delete pDimensionlessUnits;
  delete pLengthUnits;
  delete pAreaUnits;
}

void SBMLImporter::checkRuleMathConsistency ( const Rule *  pRule, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap  )

protected

Checks the expression for a give rate or assignment rule for consistency. This basically means it checks that no id present in the expression is the target for one of the following rate or assignment rules.

Definition at line 6793 of file SBMLImporter.cpp.

References CCopasiMessage::EXCEPTION, fatalError, findIdInASTTree(), getIdsFromNode(), MCSBML, mLevel, mpCopasiModel, and mVersion.

Referenced by importRuleForModelEntity().

{
  // only check if Level2 Version1
  if (this->mLevel == 2 && this->mVersion == 1)
    {
      // check if no nodes with ids of objects are used in an assignment that are
      // set in another assignment rule later on
      std::set<std::string> idSet;
      const ASTNode* pNode = pRule->getMath();
      this->getIdsFromNode(pNode, idSet);
      Model* sbmlModel = dynamic_cast<Model*>(copasi2sbmlmap[mpCopasiModel]);

      if (!sbmlModel) fatalError();

      unsigned int i, iMax = sbmlModel->getNumRules();

      for (i = 0; i < iMax; ++i)
        {
          if (sbmlModel->getRule(i) == pRule)
            {
              break;
            }
        }

      Rule* pR;
      int type;

      while (i < iMax)
        {
          pR = sbmlModel->getRule(i);
          type = pR->getTypeCode();

          if (type == SBML_ASSIGNMENT_RULE)
            {
              if (idSet.find(dynamic_cast<AssignmentRule*>(pR)->getVariable()) != idSet.end())
                {
                  CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 37, dynamic_cast<AssignmentRule*>(pR)->getVariable().c_str());
                }
            }

          ++i;
        }

      // Check if there is a reference to a reaction in the expression
      // This is not allowed for L2V1
      const ASTNode* pMath = pRule->getMath();

      if (pMath != NULL)
        {
          std::set<std::string> reactionIds;

          for (i = 0; i < sbmlModel->getListOfReactions()->size(); i++)
            {
              reactionIds.insert(sbmlModel->getReaction(i)->getId());
            }

          std::string id = SBMLImporter::findIdInASTTree(pMath, reactionIds);

          if (!id.empty())
            {
              CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 81, id.c_str());
            }
        }
    }

  // In SBML Level 3 documents,
  if (this->mLevel == 3)
    {
      // TODO we need to check if we have a species reference as the target of
      // TODO a rule because COPASI currently can not handle this and we have
      // TODO to warn the user that we will ignore the rule
      //
      // TODO Likewise we have to check if the id of a species reference is reference in the
      // TODO expression of the rule because this is also not implemented in COPASI yet and
      // TODO we ignore the rule
    }
}

bool SBMLImporter::containsVolume ( const ASTNode *  pNode, const std::string &  compartmentCN  )

protected

Checks if the volume with the given CN is one of the parameters to the given call node.

Definition at line 6065 of file SBMLImporter.cpp.

Referenced by createCReactionFromReaction().

{
  bool result = false;
  unsigned int i, iMax = pNode->getNumChildren();

  for (i = 0; i < iMax; ++i)
    {
      if (pNode->getChild(i)->getType() == AST_NAME &&
          pNode->getChild(i)->getName() == compartmentSBMLId)
        {
          result = true;
          break;
        }
    }

  return result;
}

Unit * SBMLImporter::convertSBMLCubicmetresToLitres ( const Unit *  pU )

staticprotected

If the given UnitDefinition can be converted to a form of litre, the function return the UnitDefinition in litre, otherwise NULL is returned.

Definition at line 8361 of file SBMLImporter.cpp.

References normalizeSBMLUnit(), and pResult.

Referenced by handleVolumeUnit().

{
  Unit* pResult = NULL;

  if (pU != NULL)
    {
      if ((pU->getKind() == UNIT_KIND_METER || pU->getKind() == UNIT_KIND_METRE) && pU->getExponent() % 3 == 0)
        {
          pResult = dynamic_cast<Unit*>(pU->clone());
          assert(pResult != NULL);
          Unit::removeScale(pResult);
          pResult->setExponent(pResult->getExponent() / 3);
          pResult->setKind(UNIT_KIND_LITRE);
          pResult->setMultiplier(pow(pResult->getMultiplier(), 3));
          normalizeSBMLUnit(pResult);
        }
    }

  return pResult;
}

std::map< std::string, ASTNode * > SBMLImporter::createBVarMap ( const ASTNode *  uDefFunction, const ASTNode *  function  )

protected

Creates a map of each parameter of the function definition and its corresponding parameter in the function call.

Definition at line 2953 of file SBMLImporter.cpp.

References fatalError.

{
  /* the first n-1 children, where n is the number of children, of a function definition ASTnode are the
   * arguments to the function. These correspond to the m=n-1 children of the
   * function call.
   */
  if (uDefFunction->getNumChildren() != function->getNumChildren() + 1)
    {
      std::string functionName = uDefFunction->getName();
      fatalError();
    }

  std::map<std::string, ASTNode*> varMap;
  unsigned int counter;

  for (counter = 0; counter < uDefFunction->getNumChildren() - 1; counter++)
    {
      varMap[uDefFunction->getChild(counter)->getName()] = function->getChild(counter);
    }

  return varMap;
}

CCompartment * SBMLImporter::createCCompartmentFromCompartment ( const Compartment *  sbmlCompartment, CModel *  copasiModel, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap, const Model *  pSBMLModel  )

protected

Creates and returns a COPASI CCompartment from the SBML Compartment given as argument.

Definition at line 1950 of file SBMLImporter.cpp.

References C_INVALID_INDEX, CModel::createCompartment(), CModel::getCompartments(), CCopasiVectorN< CType >::getIndex(), CCopasiContainer::getUnits(), importMIRIAM(), importNotes(), MCSBML, mLevel, CCompartment::setDimensionality(), CModelEntity::setSBMLId(), and CCopasiMessage::WARNING.

Referenced by createCModelFromSBMLDocument().

{
  if (sbmlCompartment->isSetUnits())
    {
      std::string cU = sbmlCompartment->getUnits();
    }

  unsigned int dimensionality = 3;

  if (sbmlCompartment->isSetSpatialDimensions())
    {
      dimensionality = sbmlCompartment->getSpatialDimensions();

      // starting with SBML Level 3, the spatial dimensions of a compartment can be
      // any rational number
      double dDimensionality = sbmlCompartment->getSpatialDimensions();

      if (this->mLevel > 2)
        {
          dDimensionality = sbmlCompartment->getSpatialDimensionsAsDouble();
        }

      // check if the unsigned int dimensionality corresponds to the double
      // dimensionality, otherwise give an error message
      // Actually if we wanted to be correct, we would have to check if the
      // part before the komma also matches and maybe have a separate error message if not
      dDimensionality -= dimensionality;

      if (fabs(dDimensionality) > 1e-9)
        {
          CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 89, sbmlCompartment->getId().c_str());
          dimensionality = 3;
        }
    }
  else
    {
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 90, sbmlCompartment->getId().c_str());
      dimensionality = 3;
    }

  if (dimensionality > 3)
    {
      CCopasiMessage Message(CCopasiMessage::WARNING,
                             "Reaction with id \"%s\" has dimensions of %d, this is not supported by COPASI. COPASI will assume that the compartment is three dimensional."
                             , sbmlCompartment->getId().c_str(), dimensionality);
      dimensionality = 3;
      //fatalError();
    }

  std::string name = sbmlCompartment->getName();

  if (name == "")
    {
      name = sbmlCompartment->getId();
    }

  std::string appendix = "";
  unsigned int counter = 2;
  std::ostringstream numberStream;

  while (copasiModel->getCompartments().getIndex(name + appendix) != C_INVALID_INDEX)
    {
      numberStream.str("");
      numberStream << "_" << counter;
      counter++;
      appendix = numberStream.str();
    }

  double value;
  CCompartment* copasiCompartment = copasiModel->createCompartment(name + appendix, value);

  if (this->mLevel == 1)
    {
      copasiCompartment->setSBMLId(sbmlCompartment->getName());
    }
  else
    {
      copasiCompartment->setSBMLId(sbmlCompartment->getId());
    }

  // set the dimension of the compartment
  copasiCompartment->setDimensionality(dimensionality);

  //DebugFile << "Created Compartment: " << copasiCompartment->getObjectName() << std::endl;
  SBMLImporter::importMIRIAM(sbmlCompartment, copasiCompartment);
  SBMLImporter::importNotes(copasiCompartment, sbmlCompartment);
  copasi2sbmlmap[copasiCompartment] = const_cast<Compartment*>(sbmlCompartment);
  return copasiCompartment;
}

CFunction * SBMLImporter::createCFunctionFromFunctionDefinition ( const FunctionDefinition *  sbmlFunction, CFunctionDB *  pTmpFunctionDB, Model *  pSBMLModel, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap  )

protected

Creates and returns a COPASI CFunction from the SBML FunctionDefinition given as argument.

Definition at line 1727 of file SBMLImporter.cpp.

References CFunctionDB::add(), addToKnownFunctionToMap(), areEqualFunctions(), createCFunctionFromFunctionTree(), ensureAllArgsAreBeingUsedInFunctionDefinition(), CCopasiMessage::EXCEPTION, CFunctionDB::findFunction(), functionDB, CFunction::getSBMLId(), CEvaluationTree::getType(), importMIRIAM(), importNotes(), CFunctionDB::loadedFunctions(), MCSBML, mKnownCustomUserDefinedFunctions, pdelete, CCopasiObject::setObjectName(), CFunction::setSBMLId(), CCopasiVector< T >::size(), and CEvaluationTree::UserDefined.

Referenced by importFunctionDefinitions().

{

  ensureAllArgsAreBeingUsedInFunctionDefinition(sbmlFunction);

  addToKnownFunctionToMap(mKnownCustomUserDefinedFunctions, sbmlFunction);

  CFunction* pTmpFunction = this->createCFunctionFromFunctionTree(sbmlFunction, pSBMLModel, copasi2sbmlmap);

  if (pTmpFunction == NULL)
    {
      CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 14, sbmlFunction->getId().c_str());
      return NULL;
    }

  std::string sbmlId = sbmlFunction->getId();
  pTmpFunction->setSBMLId(sbmlId);
  // check if the id is already taken by another function definition, maybe
  // from an earlier import, if this is the case, delete the id on the old
  // function definition
  // if we don't do this, two functions might have the same SBML id during
  // export which makes the exporter code so much more difficult
  size_t i, iMax = this->functionDB->loadedFunctions().size();

  for (i = 0; i < iMax; ++i)
    {
      CFunction* pFun = this->functionDB->loadedFunctions()[i];

      if (pFun->getSBMLId() == sbmlId)
        {
          pFun->setSBMLId("");
        }
    }

  std::string functionName = sbmlFunction->getName();

  if (functionName == "")
    {
      functionName = sbmlFunction->getId();
    }

  unsigned int counter = 1;
  std::ostringstream numberStream;
  std::string appendix = "";
  CFunction * pExistingFunction = NULL;

  while ((pExistingFunction = functionDB->findFunction(functionName + appendix)))
    {
      if (areEqualFunctions(pExistingFunction, pTmpFunction))
        {
          pdelete(pTmpFunction);
          pTmpFunction = pExistingFunction;

          break;
        }

      // We need to check whether the functions are identical.
      numberStream.str("");
      numberStream << "_" << counter;
      counter++;
      appendix = numberStream.str();
    }

  if (pTmpFunction != pExistingFunction)
    {
      pTmpFunction->setObjectName(functionName + appendix);
      functionDB->add(pTmpFunction, true);
      pTmpFunctionDB->add(pTmpFunction, false);
    }

  if (pTmpFunction->getType() == CEvaluationTree::UserDefined)
    {
      SBMLImporter::importMIRIAM(sbmlFunction, pTmpFunction);
      SBMLImporter::importNotes(pTmpFunction, sbmlFunction);
    }

  return pTmpFunction;
}

CFunction * SBMLImporter::createCFunctionFromFunctionTree ( const FunctionDefinition *  pSBMLFunction, Model *  pSBMLModel, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap  )

protected

Definition at line 1806 of file SBMLImporter.cpp.

References CFunction::addVariable(), CEvaluationNodeVariable::ANY, CEvaluationTree::compile(), CCopasiMessage::EXCEPTION, CEvaluationTree::getRoot(), isDelayFunctionUsed(), MCSBML, mExplicitelyTimeDependentFunctionDefinitions, preprocessNode(), replaceTimeNodesInFunctionDefinition(), CEvaluationTree::setRoot(), CEvaluationTree::setTree(), and CEvaluationTree::updateTree().

Referenced by createCFunctionFromFunctionDefinition().

{
  CFunction* pFun = NULL;

  if (!pSBMLFunction->isSetMath())
    return NULL;

  ConverterASTNode root(*pSBMLFunction->getMath());

  if (SBMLImporter::isDelayFunctionUsed(&root))
    {
      CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 85, pSBMLFunction->getId().c_str());
    }

  this->preprocessNode(&root, pSBMLModel, copasi2sbmlmap);

  if (root.getType() != AST_LAMBDA)
    {
      CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 11, pSBMLFunction->getId().c_str());
      return NULL;
    }

  // get the number of children.
  // the first n-1 children are the parameters for the function
  // the last child is the actual function
  pFun = new CKinFunction();
  unsigned int i, iMax = root.getNumChildren() - 1;
  std::set<std::string> variableNames;

  for (i = 0; i < iMax; ++i)
    {
      ASTNode* pVarNode = root.getChild(i);

      if (pVarNode->getType() != AST_NAME)
        {
          delete pFun;
          CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 12, pSBMLFunction->getId().c_str());
        }

      pFun->addVariable(pVarNode->getName());
      variableNames.insert(pVarNode->getName());
    }

  // now we check if the AST tree has a node that represents the time
  // object
  // find a unique name for the time variable
  std::ostringstream sstream;
  std::string timeVariableName = "time";
  unsigned int postfix = 1;

  while (variableNames.find(timeVariableName) != variableNames.end())
    {
      sstream.str("");
      sstream << "time_" << postfix;
      timeVariableName = sstream.str();
      ++postfix;
    }

  if (this->replaceTimeNodesInFunctionDefinition(root.getChild(iMax), timeVariableName))
    {
      // add another variable to the function
      ASTNode* pVarNode = new ASTNode(AST_NAME);
      pVarNode->setName(timeVariableName.c_str());
      ASTNode* pTmpNode = root.removeChild(iMax);
      root.addChild(pVarNode);
      root.addChild(pTmpNode);
      // increase iMax since we now have one more child
      ++iMax;
      pFun->addVariable(timeVariableName);
      this->mExplicitelyTimeDependentFunctionDefinitions.insert(pSBMLFunction->getId());
    }

  pFun->setTree(*root.getChild(iMax));
  CCopasiTree<CEvaluationNode>::iterator treeIt = pFun->getRoot();

  // if the root node already is an object node, this has to be dealt with separately
  if (dynamic_cast<CEvaluationNodeObject*>(&(*treeIt)))
    {
      CEvaluationNodeVariable* pVariableNode = new CEvaluationNodeVariable(CEvaluationNodeVariable::ANY, (*treeIt).getData().substr(1, (*treeIt).getData().length() - 2));
      pFun->setRoot(pVariableNode);
    }
  else
    {
      while (treeIt != NULL)
        {
          if (dynamic_cast<CEvaluationNodeObject*>(&(*treeIt)))
            {
              CEvaluationNodeVariable* pVariableNode = new CEvaluationNodeVariable(CEvaluationNodeVariable::ANY, (*treeIt).getData().substr(1, (*treeIt).getData().length() - 2));

              if ((*treeIt).getParent())
                {
                  (*treeIt).getParent()->addChild(pVariableNode, &(*treeIt));
                  (*treeIt).getParent()->removeChild(&(*treeIt));
                }

              delete &(*treeIt);
              treeIt = pVariableNode;
            }

          ++treeIt;
        }
    }

  pFun->updateTree();

  if (!pFun->compile())
    {
      delete pFun;
      CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 28, pSBMLFunction->getId().c_str());
    }

  if (pFun->getRoot() == NULL)
    {
      delete pFun;
      CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 13, pSBMLFunction->getId().c_str());
    }

  return pFun;
}

CMetab * SBMLImporter::createCMetabFromSpecies ( const Species *  sbmlSpecies, CModel *  copasiModel, CCompartment *  copasiCompartment, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap, const Model *  pSBMLModel  )

protected

Creates and returns a COPASI CMetab from the given SBML Species object.

Creates and returns a Copasi CMetab from the given SBML Species object.

Definition at line 2044 of file SBMLImporter.cpp.

References C_INVALID_INDEX, CModel::createMetabolite(), fatalError, CModelEntity::FIXED, CCopasiVectorN< CType >::getIndex(), CCompartment::getMetabolites(), CCopasiObject::getObjectName(), if(), importMIRIAM(), importNotes(), mLevel, mSubstanceOnlySpecies, CModelEntity::REACTIONS, CModelEntity::setSBMLId(), and CMetab::setStatus().

Referenced by createCModelFromSBMLDocument().

{
  if (sbmlSpecies->isSetSubstanceUnits())
    {
      std::string cU = sbmlSpecies->getSubstanceUnits();
    }

  std::map<CCopasiObject*, SBase*>::iterator it = copasi2sbmlmap.find(copasiCompartment);

  if (it == copasi2sbmlmap.end())
    {
      fatalError();
    }

  Compartment* pSBMLCompartment = (Compartment*)it->second;

  if (sbmlSpecies->getHasOnlySubstanceUnits() == true)
    {
      this->mSubstanceOnlySpecies.insert(std::make_pair(const_cast<Species*>(sbmlSpecies), pSBMLCompartment));
    }

  std::string name = sbmlSpecies->getName();

  if (name == "")
    {
      name = sbmlSpecies->getId();
    }

  std::string appendix = "";
  unsigned int counter = 2;
  std::ostringstream numberStream;

  while (copasiCompartment->getMetabolites().getIndex(name + appendix) != C_INVALID_INDEX)
    {
      numberStream.str("");
      numberStream << "_" << counter;
      counter++;
      appendix = numberStream.str();
    }

  CMetab* copasiMetabolite = copasiModel->createMetabolite(name + appendix, copasiCompartment->getObjectName());

  if (copasiMetabolite == NULL)
    {
      //DebugFile << "Could not create Copasi species." << std::endl;
      fatalError();
    }

  if (sbmlSpecies->getConstant() || sbmlSpecies->getBoundaryCondition())
    {
      copasiMetabolite->setStatus(CModelEntity::FIXED);
    }
  else
    {
      copasiMetabolite->setStatus(CModelEntity::REACTIONS);
    }

  // also check if the compartment has a spatialSize of 0 because this also implies hasOnlySubstanceUnits for the species in this compartment
  // In Level 3 files this no longer seems to be the case however, so we only
  // do this for L1 and L2 files
  if (pSBMLCompartment->isSetSpatialDimensions() && pSBMLCompartment->getSpatialDimensions() == 0 && this->mLevel < 3)
    {
      this->mSubstanceOnlySpecies.insert(std::make_pair(const_cast<Species*>(sbmlSpecies), pSBMLCompartment));
    }

  //DebugFile << "Created species: " << copasiMetabolite->getObjectName() << std::endl;
  copasi2sbmlmap[copasiMetabolite] = const_cast<Species*>(sbmlSpecies);

  if (this->mLevel == 1)
    {
      copasiMetabolite->setSBMLId(sbmlSpecies->getName());
    }
  else
    {
      copasiMetabolite->setSBMLId(sbmlSpecies->getId());
    }

  SBMLImporter::importMIRIAM(sbmlSpecies, copasiMetabolite);
  SBMLImporter::importNotes(copasiMetabolite, sbmlSpecies);
  return copasiMetabolite;
}

CModel * SBMLImporter::createCModelFromSBMLDocument ( SBMLDocument *  sbmlDocument, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap  )

protected

Creates and returns a COPASI CModel from the SBMLDocument given as argument.

Creates and returns a Copasi CModel from the SBMLDocument given as argument.

Definition at line 243 of file SBMLImporter.cpp.

References CProcessReport::addItem(), applyStoichiometricExpressions(), areRulesUnique(), areSBMLUnitDefinitionsIdentical(), C_INT32, C_INVALID_INDEX, checkElementUnits(), createCCompartmentFromCompartment(), createCMetabFromSpecies(), createCModelValueFromParameter(), createCReactionFromReaction(), CModel::deterministic, CCopasiMessage::EXCEPTION, fatalError, findAvogadroConstant(), CProcessReport::finishItem(), CModel::forceCompile(), functionDB, getInitialCNForSBase(), CCopasiContainer::getObject(), getOriginalSBMLId(), CModel::getQuantity2NumberFactor(), CFunction::getSBMLId(), CCopasiMessage::getText(), handleAreaUnit(), handleLengthUnit(), handleSubstanceUnit(), handleTimeUnit(), handleVolumeUnit(), importEvents(), importFunctionDefinitions(), importInitialAssignments(), importMIRIAM(), importNotes(), importSBMLRule(), INIT_DEFAULTS, isStochasticModel(), CModel::l, CFunctionDB::loadedFunctions(), CModel::m, CModel::m2, mAvogadroCreated, mAvogadroSet, MCSBML, mDelayFound, mDivisionByCompartmentReactions, mFastReactions, mhImportStep, mIgnoredParameterUnits, mImportStep, mKnownInitalValues, mLevel, CModel::Mol, mpCopasiModel, mpDataModel, mpImportHandler, mPotentialAvogadroNumbers, mReactionsWithReplacedLocalParameters, mSubstanceOnlySpecies, mUnitOnNumberFound, mUnsupportedRuleFound, mVersion, CCopasiMessage::peekLastMessage(), CProcessReport::progressItem(), removeUnusedFunctions(), CModel::s, sbmlIdMap, CModel::setAreaUnit(), CModel::setAvogadro(), CModel::setCompileFlag(), setInitialValues(), CModel::setLengthUnit(), CModel::setModelType(), CCopasiObject::setObjectName(), CModel::setQuantityUnit(), CFunction::setSBMLId(), CModelEntity::setSBMLId(), CModel::setTimeUnit(), CModel::setVolumeUnit(), CCopasiVector< T >::size(), speciesMap, CModel::stochastic, and CCopasiMessage::WARNING.

Referenced by parseSBML().

{
  Model* sbmlModel = sbmlDocument->getModel();

  /* Create an empty model and set the title. */
  this->mpCopasiModel = new CModel(mpDataModel);
  copasi2sbmlmap[this->mpCopasiModel] = sbmlModel;
  this->mpCopasiModel->setLengthUnit(CModel::m);
  this->mpCopasiModel->setAreaUnit(CModel::m2);
  this->mpCopasiModel->setVolumeUnit(CModel::l);
  this->mpCopasiModel->setTimeUnit(CModel::s);
  this->mpCopasiModel->setQuantityUnit(CModel::Mol);
  this->mpCopasiModel->setSBMLId(sbmlModel->getId());

  unsigned C_INT32 step = 0, totalSteps = 0;
  size_t hStep = C_INVALID_INDEX;

  mImportStep = 1;

  if (mpImportHandler && !mpImportHandler->progressItem(mhImportStep)) return NULL;

  SBMLImporter::importMIRIAM(sbmlModel, this->mpCopasiModel);
  UnitDefinition *pSubstanceUnits = NULL;
  UnitDefinition *pTimeUnits = NULL;
  UnitDefinition *pVolumeUnits = NULL;
  UnitDefinition *pAreaUnits = NULL;
  UnitDefinition *pLengthUnits = NULL;
  /* Set standard units to match the standard units of SBML files. */

  // for SBML L3 files the default units are defined on the model
  if (this->mLevel > 2)
    {
      this->mpCopasiModel->setAvogadro(6.02214179e23);
      this->mAvogadroSet = true;

      // we make copies of the unit definitions so that we do not have to remember
      // if we created them or not
      std::string units;
      UnitDefinition* pUDef = NULL;
      Unit unit(sbmlModel->getLevel(), sbmlModel->getVersion());
      INIT_DEFAULTS(unit);

      if (sbmlModel->isSetSubstanceUnits())
        {
          units = sbmlModel->getSubstanceUnits();
          assert(units != "");
          pUDef =  sbmlModel->getUnitDefinition(units);

          if (pUDef != NULL)
            {
              pSubstanceUnits = new UnitDefinition(*pUDef);
            }
          else
            {
              if (units == "mole")
                {
                  unit.setKind(UNIT_KIND_MOLE);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pSubstanceUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pSubstanceUnits->addUnit(&unit);
                }
              else if (units == "item")
                {
                  unit.setKind(UNIT_KIND_ITEM);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pSubstanceUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pSubstanceUnits->addUnit(&unit);
                }
              else if (units == "dimensionless")
                {
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pSubstanceUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pSubstanceUnits->addUnit(&unit);
                }
              else
                {
                  std::string message = "COPASI can't handle substance unit \"" + units + "\". Setting unit for substances to dimensionless.";
                  CCopasiMessage(CCopasiMessage::WARNING, message.c_str());
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pSubstanceUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pSubstanceUnits->addUnit(&unit);
                }
            }
        }

      if (sbmlModel->isSetTimeUnits())
        {
          units = sbmlModel->getTimeUnits();
          assert(units != "");
          pUDef =  sbmlModel->getUnitDefinition(units);

          if (pUDef != NULL)
            {
              pTimeUnits = new UnitDefinition(*pUDef);
            }
          else
            {
              if (units == "second")
                {
                  unit.setKind(UNIT_KIND_SECOND);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pTimeUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pTimeUnits->addUnit(&unit);
                }
              else if (units == "dimensionless")
                {
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pTimeUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pTimeUnits->addUnit(&unit);
                }
              else
                {
                  std::string message = "COPASI can't handle time unit \"" + units + "\". Setting unit for time to dimensionless.";
                  CCopasiMessage(CCopasiMessage::WARNING, message.c_str());
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pTimeUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pTimeUnits->addUnit(&unit);
                }
            }
        }

      if (sbmlModel->isSetVolumeUnits())
        {
          units = sbmlModel->getVolumeUnits();
          assert(units != "");
          pUDef =  sbmlModel->getUnitDefinition(units);

          if (pUDef != NULL)
            {
              pVolumeUnits = new UnitDefinition(*pUDef);
            }
          else
            {
              if (units == "litre")
                {
                  unit.setKind(UNIT_KIND_LITRE);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pVolumeUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pVolumeUnits->addUnit(&unit);
                }
              else if (units == "dimensionless")
                {
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pVolumeUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pVolumeUnits->addUnit(&unit);
                }
              else
                {
                  std::string message = "COPASI can't handle volume unit \"" + units + "\". Setting unit for volume to dimensionless.";
                  CCopasiMessage(CCopasiMessage::WARNING, message.c_str());
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pVolumeUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pVolumeUnits->addUnit(&unit);
                }
            }
        }

      if (sbmlModel->isSetAreaUnits())
        {
          units = sbmlModel->getAreaUnits();
          assert(units != "");
          pUDef =  sbmlModel->getUnitDefinition(units);

          if (pUDef != NULL)
            {
              pAreaUnits = new UnitDefinition(*pUDef);
            }
          else
            {
              if (units == "dimensionless")
                {
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pAreaUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pAreaUnits->addUnit(&unit);
                }
              else
                {
                  std::string message = "COPASI can't handle area unit \"" + units + "\". Setting unit for area to dimensionless.";
                  CCopasiMessage(CCopasiMessage::WARNING, message.c_str());
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pAreaUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pAreaUnits->addUnit(&unit);
                }
            }
        }

      if (sbmlModel->isSetLengthUnits())
        {
          units = sbmlModel->getLengthUnits();
          assert(units != "");
          pUDef =  sbmlModel->getUnitDefinition(units);

          if (pUDef != NULL)
            {
              pLengthUnits = new UnitDefinition(*pUDef);
            }
          else
            {
              if (units == "litre")
                {
                  unit.setKind(UNIT_KIND_LITRE);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pLengthUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pLengthUnits->addUnit(&unit);
                }
              else if (units == "metre")
                {
                  unit.setKind(UNIT_KIND_METRE);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pLengthUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pLengthUnits->addUnit(&unit);
                }
              else if (units == "dimensionless")
                {
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pLengthUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pLengthUnits->addUnit(&unit);
                }
              else
                {
                  std::string message = "COPASI can't handle length unit \"" + units + "\". Setting unit for length to dimensionless.";
                  CCopasiMessage(CCopasiMessage::WARNING, message.c_str());
                  unit.setKind(UNIT_KIND_DIMENSIONLESS);
                  unit.setExponent(1);
                  unit.setMultiplier(1.0);
                  unit.setScale(0);
                  pLengthUnits = new UnitDefinition(sbmlModel->getLevel(), sbmlModel->getVersion());
                  pLengthUnits->addUnit(&unit);
                }
            }
        }
    }
  else
    {
      if (sbmlModel->getNumUnitDefinitions() != 0)
        {
          unsigned int counter;

          // we make copies o the unit definitions so that we do not have to remember
          // if we created them or not
          for (counter = 0; counter < sbmlModel->getNumUnitDefinitions(); counter++)
            {
              UnitDefinition* uDef = sbmlModel->getUnitDefinition(counter);
              std::string unitId = uDef->getId();

              if (unitId == "substance")
                {
                  pSubstanceUnits = new UnitDefinition(*uDef);
                }
              else if (unitId == "time")
                {
                  pTimeUnits = new UnitDefinition(*uDef);
                }
              else if (unitId == "volume")
                {
                  pVolumeUnits = new UnitDefinition(*uDef);
                }
              else if ((unitId == "area"))
                {
                  pAreaUnits = new UnitDefinition(*uDef);
                }
              else if ((unitId == "length"))
                {
                  pLengthUnits = new UnitDefinition(*uDef);
                }
            }
        }
    }

  // create the default units if some unit has not been specified
  if (pSubstanceUnits == NULL)
    {
      if (this->mLevel > 2)
        {
          // issue a warning
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 91, "substance", "mole" , "substance");
        }

      // create the default units
      pSubstanceUnits = new UnitDefinition(this->mLevel, this->mVersion);
      pSubstanceUnits->setId("dummy_substance");
      pSubstanceUnits->setName("dummy_substance");
      Unit* pUnit = pSubstanceUnits->createUnit();
      pUnit->setKind(UNIT_KIND_MOLE);
      pUnit->setExponent(1);
      pUnit->setMultiplier(1.0);
      pUnit->setScale(0);
    }

  if (pTimeUnits == NULL)
    {
      if (this->mLevel > 2)
        {
          // issue a warning
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 91, "time", "second" , "time");
        }

      // create the default units
      pTimeUnits = new UnitDefinition(this->mLevel, this->mVersion);
      pTimeUnits->setId("dummy_time");
      pTimeUnits->setName("dummy_time");
      Unit* pUnit = pTimeUnits->createUnit();
      pUnit->setKind(UNIT_KIND_SECOND);
      pUnit->setExponent(1);
      pUnit->setMultiplier(1.0);
      pUnit->setScale(0);
    }

  if (pVolumeUnits == NULL)
    {
      if (this->mLevel > 2)
        {
          // issue a warning
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 91, "volume", "litre" , "volume");
        }

      // create the default units
      pVolumeUnits = new UnitDefinition(this->mLevel, this->mVersion);
      pVolumeUnits->setId("dummy_volume");
      pVolumeUnits->setName("dummy_volume");
      Unit* pUnit = pVolumeUnits->createUnit();
      pUnit->setKind(UNIT_KIND_LITRE);
      pUnit->setExponent(1);
      pUnit->setMultiplier(1.0);
      pUnit->setScale(0);
    }

  if (pAreaUnits == NULL)
    {
      if (this->mLevel > 2)
        {
          // issue a warning
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 91, "area", "m^2" , "area");
        }

      // create the default units
      pAreaUnits = new UnitDefinition(this->mLevel, this->mVersion);
      pAreaUnits->setId("dummy_area");
      pAreaUnits->setName("dummy_area");
      Unit* pUnit = pAreaUnits->createUnit();
      pUnit->setKind(UNIT_KIND_METRE);
      pUnit->setExponent(2);
      pUnit->setMultiplier(1.0);
      pUnit->setScale(0);
    }

  if (pLengthUnits == NULL)
    {
      if (this->mLevel > 2)
        {
          // issue a warning
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 91, "length", "m" , "length");
        }

      // create the default units
      pLengthUnits = new UnitDefinition(this->mLevel, this->mVersion);
      pLengthUnits->setId("dummy_length");
      pLengthUnits->setName("dummy_length");
      Unit* pUnit = pLengthUnits->createUnit();
      pUnit->setKind(UNIT_KIND_METRE);
      pUnit->setExponent(1);
      pUnit->setMultiplier(1.0);
      pUnit->setScale(0);
    }

  // now we have some common code to actually import the units

  // handle the substance units
  assert(pSubstanceUnits != NULL);

  if (pSubstanceUnits != NULL)
    {
      std::pair<CModel::QuantityUnit, bool> qUnit;

      try
        {
          qUnit = this->handleSubstanceUnit(pSubstanceUnits);
        }
      catch (...)
        {
          std::ostringstream os;
          os << "Error while importing substance units.";

          // check if the last message on the stack is an exception
          // and if so, add the message text to the current exception
          if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
            {
              // we only want the message, not the timestamp line
              std::string text = CCopasiMessage::peekLastMessage().getText();
              os << text.substr(text.find("\n"));
            }

          CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
        }

      if (qUnit.second == false)
        {
          // the unit could not be handled, give an error message and
          // set the units to mole
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 66, "substance", "Mole");
          this->mpCopasiModel->setQuantityUnit(CModel::Mol);
        }
      else
        {
          this->mpCopasiModel->setQuantityUnit(qUnit.first);
        }

      // check if the extends units are set and if they are equal to the substance units
      // otherwise issue a warning
      if (this->mLevel > 2)
        {
          if (!sbmlModel->isSetExtentUnits())
            {
              CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 92);
            }
          else
            {
              const UnitDefinition* pExtendsUnits = sbmlModel->getUnitDefinition(sbmlModel->getExtentUnits());

              if (pExtendsUnits != NULL)
                {
                  if (!areSBMLUnitDefinitionsIdentical(pSubstanceUnits, pExtendsUnits))
                    {
                      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 92);
                    }
                }
              else
                {
                  CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 66, "extends", "the same units as the substances");
                }
            }
        }

      delete pSubstanceUnits;
      pSubstanceUnits = NULL;
    }

  // handle the time units
  assert(pTimeUnits != NULL);

  if (pTimeUnits != NULL)
    {
      std::pair<CModel::TimeUnit, bool> tUnit;

      try
        {
          tUnit = this->handleTimeUnit(pTimeUnits);
        }
      catch (...)
        {
          std::ostringstream os;
          os << "Error while importing time units.";

          // check if the last message on the stack is an exception
          // and if so, add the message text to the current exception
          if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
            {
              // we only want the message, not the timestamp line
              std::string text = CCopasiMessage::peekLastMessage().getText();
              os << text.substr(text.find("\n"));
            }

          CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
        }

      if (tUnit.second == false)
        {
          // the unit could not be handled, give an error message and
          // set the units to second
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 66, "time", "second");
          this->mpCopasiModel->setTimeUnit(CModel::s);
        }
      else
        {
          this->mpCopasiModel->setTimeUnit(tUnit.first);
        }

      delete pTimeUnits;
      pTimeUnits = NULL;
    }

  // handle the volume units
  assert(pVolumeUnits != NULL);

  if (pVolumeUnits != NULL)
    {
      std::pair<CModel::VolumeUnit, bool> vUnit;

      try
        {
          vUnit = this->handleVolumeUnit(pVolumeUnits);
        }
      catch (...)
        {
          std::ostringstream os;
          os << "Error while importing volume units.";

          // check if the last message on the stack is an exception
          // and if so, add the message text to the current exception
          if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
            {
              // we only want the message, not the timestamp line
              std::string text = CCopasiMessage::peekLastMessage().getText();
              os << text.substr(text.find("\n"));
            }

          CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
        }

      if (vUnit.second == false)
        {
          // the unit could not be handled, give an error message and
          // set the units to litre
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 66, "volume", "litre");
          this->mpCopasiModel->setVolumeUnit(CModel::l);
        }
      else
        {
          this->mpCopasiModel->setVolumeUnit(vUnit.first);
        }

      delete pVolumeUnits;
      pVolumeUnits = NULL;
    }

  // handle the area units
  assert(pAreaUnits != NULL);

  if (pAreaUnits != NULL)
    {
      std::pair<CModel::AreaUnit, bool> vUnit;

      try
        {
          vUnit = this->handleAreaUnit(pAreaUnits);
        }
      catch (...)
        {
          std::ostringstream os;
          os << "Error while importing area units.";

          // check if the last message on the stack is an exception
          // and if so, add the message text to the current exception
          if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
            {
              // we only want the message, not the timestamp line
              std::string text = CCopasiMessage::peekLastMessage().getText();
              os << text.substr(text.find("\n"));
            }

          CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
        }

      if (vUnit.second == false)
        {
          // the unit could not be handled, give an error message and
          // set the units to litre
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 66, "area", "square meter");
          this->mpCopasiModel->setAreaUnit(CModel::m2);
        }
      else
        {
          this->mpCopasiModel->setAreaUnit(vUnit.first);
        }

      delete pAreaUnits;
      pAreaUnits = NULL;
    }

  // handle the length units
  assert(pLengthUnits != NULL);

  if (pLengthUnits != NULL)
    {
      std::pair<CModel::LengthUnit, bool> vUnit;

      try
        {
          vUnit = this->handleLengthUnit(pLengthUnits);
        }
      catch (...)
        {
          std::ostringstream os;
          os << "Error while importing length units.";

          // check if the last message on the stack is an exception
          // and if so, add the message text to the current exception
          if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
            {
              // we only want the message, not the timestamp line
              std::string text = CCopasiMessage::peekLastMessage().getText();
              os << text.substr(text.find("\n"));
            }

          CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
        }

      if (vUnit.second == false)
        {
          // the unit could not be handled, give an error message and
          // set the units to litre
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 66, "length", "meter");
          this->mpCopasiModel->setLengthUnit(CModel::m);
        }
      else
        {
          this->mpCopasiModel->setLengthUnit(vUnit.first);
        }

      delete pLengthUnits;
      pLengthUnits = NULL;
    }

  // go through all compartments and species and check if the units are
  // consistent
  checkElementUnits(sbmlModel, this->mpCopasiModel, this->mLevel, this->mVersion);
  std::string title;

  if (this->isStochasticModel(sbmlModel))
    {
      this->mpCopasiModel->setModelType(CModel::stochastic);
    }
  else
    {
      this->mpCopasiModel->setModelType(CModel::deterministic);
    }

  SBMLImporter::importNotes(this->mpCopasiModel, sbmlModel);

  title = sbmlModel->getName();

  if (title == "")
    {
      title = "NoName";
    }

  size_t idCount = 0;

  while (mpDataModel->getObject("Model=" + title) != NULL)
    {
      std::stringstream str; str << sbmlModel->getName() << "_" << ++idCount;
      title = str.str();
    }

  this->mpCopasiModel->setObjectName(title.c_str());
  // fill the set of SBML species reference ids because
  // we need this to check for references to species references in all expressions
  // as long as we do not support these references
  SBMLImporter::updateSBMLSpeciesReferenceIds(sbmlModel, this->mSBMLSpeciesReferenceIds);

  /* import the functions */
  unsigned int counter;
  CCopasiVectorN< CFunction >* functions = &(this->functionDB->loadedFunctions());

  size_t num = (*functions).size();

  if (mpImportHandler)
    {
      mImportStep = 2;

      if (!mpImportHandler->progressItem(mhImportStep)) return NULL;

      totalSteps = (unsigned C_INT32) num;
      hStep = mpImportHandler->addItem("Importing function definitions",
                                       step,
                                       &totalSteps);
    }

  this->sbmlIdMap.clear();

  for (counter = 0; counter < num; ++counter)
    {
      CFunction* tree = (*functions)[counter];

      if (!tree->getSBMLId().empty())
        {
          this->sbmlIdMap[tree] = tree->getSBMLId();
          tree->setSBMLId("");
        }
    }

  CFunctionDB* pTmpFunctionDB = this->importFunctionDefinitions(sbmlModel, copasi2sbmlmap);
  // try to find global parameters that represent avogadros number
  this->findAvogadroConstant(sbmlModel, this->mpCopasiModel->getQuantity2NumberFactor());

  std::map<std::string, CCompartment*> compartmentMap;

  /* Create the compartments */
  num = sbmlModel->getNumCompartments();

  if (mpImportHandler)
    {
      mpImportHandler->finishItem(hStep);
      mImportStep = 3;

      if (!mpImportHandler->progressItem(mhImportStep)) return NULL;

      step = 0;
      totalSteps = (unsigned C_INT32) num;
      hStep = mpImportHandler->addItem("Importing compartments...",
                                       step,
                                       &totalSteps);
    }

  for (counter = 0; counter < num; counter++)
    {
      Compartment* sbmlCompartment = sbmlModel->getCompartment(counter);

      if (sbmlCompartment == NULL)
        {
          fatalError();
        }

      CCompartment* pCopasiCompartment = NULL;

      try
        {
          pCopasiCompartment = this->createCCompartmentFromCompartment(sbmlCompartment, this->mpCopasiModel, copasi2sbmlmap, sbmlModel);
        }
      catch (...)
        {
          std::ostringstream os;
          os << "Error while importing compartment \"";
          os << sbmlCompartment->getId() << "\".";

          // check if the last message on the stack is an exception
          // and if so, add the message text to the current exception
          if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
            {
              // we only want the message, not the timestamp line
              std::string text = CCopasiMessage::peekLastMessage().getText();
              os << text.substr(text.find("\n"));
            }

          CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
        }

      std::string key = sbmlCompartment->getId();

      if (mLevel == 1)
        {
          key = sbmlCompartment->getName();
        }

      compartmentMap[key] = pCopasiCompartment;
      ++step;

      if (mpImportHandler && !mpImportHandler->progressItem(hStep)) return NULL;
    }

  /* Create all species */
  num = sbmlModel->getNumSpecies();

  if (mpImportHandler)
    {
      mpImportHandler->finishItem(hStep);
      mImportStep = 4;

      if (!mpImportHandler->progressItem(mhImportStep)) return NULL;

      step = 0;
      totalSteps = (unsigned C_INT32) num;
      hStep = mpImportHandler->addItem("Importing species...",
                                       step,
                                       &totalSteps);
    }

  for (counter = 0; counter < num; ++counter)
    {
      Species* sbmlSpecies = sbmlModel->getSpecies(counter);

      if (sbmlSpecies == NULL)
        {
          fatalError();
        }

      CCompartment* pCopasiCompartment = compartmentMap[sbmlSpecies->getCompartment()];

      if (pCopasiCompartment != NULL)
        {
          CMetab* pCopasiMetabolite = NULL;

          try
            {
              pCopasiMetabolite = this->createCMetabFromSpecies(sbmlSpecies, this->mpCopasiModel, pCopasiCompartment, copasi2sbmlmap, sbmlModel);
            }
          catch (...)
            {
              std::ostringstream os;
              os << "Error while importing species \"";
              os << sbmlSpecies->getId() << "\".";

              // check if the last message on the stack is an exception
              // and if so, add the message text to the current exception
              if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
                {
                  // we only want the message, not the timestamp line
                  std::string text = CCopasiMessage::peekLastMessage().getText();
                  os << text.substr(text.find("\n"));
                }

              CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
            }

          std::string key;

          if (this->mLevel == 1)
            {
              key = sbmlSpecies->getName();
            }
          else
            {
              key = sbmlSpecies->getId();
            }

          this->speciesMap[key] = pCopasiMetabolite;
        }
      else
        {
          CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 5 , sbmlSpecies->getCompartment().c_str(), sbmlSpecies->getId().c_str());
        }

      ++step;

      if (mpImportHandler && !mpImportHandler->progressItem(hStep)) return NULL;
    }

  /* Create the global Parameters */
  num = sbmlModel->getNumParameters();

  if (mpImportHandler)
    {
      mpImportHandler->finishItem(hStep);
      mImportStep = 5;

      if (!mpImportHandler->progressItem(mhImportStep)) return NULL;

      step = 0;
      totalSteps = (unsigned C_INT32) num;
      hStep = mpImportHandler->addItem("Importing global parameters...",
                                       step,
                                       &totalSteps);
    }

  for (counter = 0; counter < num; ++counter)
    {
      Parameter* sbmlParameter = sbmlModel->getParameter(counter);

      if (sbmlParameter == NULL)
        {
          fatalError();
        }

      try
        {
          std::string sbmlId = getOriginalSBMLId(sbmlParameter);

          // don't import parameter if it is one of our special ones instead get the cn to the target
          if (!sbmlId.empty())
            {
              Species* species = sbmlModel->getSpecies(sbmlId);

              if (species != NULL)
                {
                  mKnownInitalValues[sbmlParameter->getId()] = getInitialCNForSBase(species, copasi2sbmlmap);
                  continue;
                }

              Compartment *comp = sbmlModel->getCompartment(sbmlId);

              if (comp != NULL)
                {
                  mKnownInitalValues[sbmlParameter->getId()] = getInitialCNForSBase(comp, copasi2sbmlmap);
                  continue;
                }

              Parameter* param = sbmlModel->getParameter(sbmlId);

              if (param != NULL)
                {
                  mKnownInitalValues[sbmlParameter->getId()] = getInitialCNForSBase(param, copasi2sbmlmap);
                  continue;
                }
            }

          this->createCModelValueFromParameter(sbmlParameter, this->mpCopasiModel, copasi2sbmlmap);
        }
      catch (...)
        {
          std::ostringstream os;
          os << "Error while importing parameter \"";
          os << sbmlParameter->getId() << "\".";

          // check if the last message on the stack is an exception
          // and if so, add the message text to the current exception
          if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
            {
              // we only want the message, not the timestamp line
              std::string text = CCopasiMessage::peekLastMessage().getText();
              os << text.substr(text.find("\n"));
            }

          CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
        }

      ++step;

      if (mpImportHandler && !mpImportHandler->progressItem(hStep)) return NULL;
    }

  // the information that is collected here is used in the creation of the reactions to
  // adjust the stoichiometry of substrates and products with the conversion factos from
  // the SBML model

  // now that all parameters have been imported, we can check if the model
  // defines a global conversion factor
  if (this->mLevel > 2 && sbmlModel->isSetConversionFactor())
    {
      this->mConversionFactorFound = true;
      std::string id = sbmlModel->getConversionFactor();
      assert(id != "");
      std::map<std::string, const CModelValue*>::const_iterator pos = this->mSBMLIdModelValueMap.find(id);
      assert(pos != this->mSBMLIdModelValueMap.end());

      if (pos != this->mSBMLIdModelValueMap.end())
        {
          this->mpModelConversionFactor = pos->second;
        }
      else
        {
          fatalError();
        }
    }

  // we need to go through all species again and find conversion factors
  // right now I don't want to change the order of importing model values
  // and species since I am not sure what implications that would have
  if (this->mLevel > 2)
    {
      num = sbmlModel->getNumSpecies();

      for (counter = 0; counter < num; ++counter)
        {
          Species* sbmlSpecies = sbmlModel->getSpecies(counter);

          if (sbmlSpecies->isSetConversionFactor())
            {
              this->mConversionFactorFound = true;
              std::map<std::string, const CModelValue*>::const_iterator pos = this->mSBMLIdModelValueMap.find(sbmlSpecies->getConversionFactor());

              if (pos != this->mSBMLIdModelValueMap.end())
                {
                  this->mSpeciesConversionParameterMap[sbmlSpecies->getId()] = pos->second;
                }
            }
        }
    }

  if (!this->mIgnoredParameterUnits.empty())
    {
      std::ostringstream os;
      std::vector<std::string>::iterator errorIt = this->mIgnoredParameterUnits.begin();

      while (errorIt != this->mIgnoredParameterUnits.end())
        {
          os << *errorIt << ", ";
          ++errorIt;
        }

      std::string s = os.str();
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 26, s.substr(0, s.size() - 2).c_str());
    }

  /* Create all reactions */
  num = sbmlModel->getNumReactions();

  if (mpImportHandler)
    {
      mpImportHandler->finishItem(hStep);
      mImportStep = 6;

      if (!mpImportHandler->progressItem(mhImportStep)) return NULL;

      step = 0;
      totalSteps = (unsigned C_INT32) num;
      hStep = mpImportHandler->addItem("Importing reactions...",
                                       step,
                                       &totalSteps);
    }

  this->mDivisionByCompartmentReactions.clear();
  this->mFastReactions.clear();
  this->mReactionsWithReplacedLocalParameters.clear();

  for (counter = 0; counter < num; counter++)
    {
      try
        {
          this->createCReactionFromReaction(sbmlModel->getReaction(counter), sbmlModel, this->mpCopasiModel, copasi2sbmlmap, pTmpFunctionDB);
        }
      catch (...)
        {
          std::ostringstream os;
          os << "Error while importing reaction \"";
          os << sbmlModel->getReaction(counter)->getId() << "\".";

          // check if the last message on the stack is an exception
          // and if so, add the message text to the current exception
          if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
            {
              // we only want the message, not the timestamp line
              std::string text = CCopasiMessage::peekLastMessage().getText();
              os << text.substr(text.find("\n"));
            }

          CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
        }

      ++step;

      if (mpImportHandler && !mpImportHandler->progressItem(hStep)) return NULL;
    }

  if (!this->mDivisionByCompartmentReactions.empty())
    {
      // create the error message
      std::string idList;
      std::set<std::string>::const_iterator it = this->mDivisionByCompartmentReactions.begin();
      std::set<std::string>::const_iterator endit = this->mDivisionByCompartmentReactions.end();

      while (it != endit)
        {
          idList += (*it);
          idList += ", ";
          ++it;
        }

      idList = idList.substr(0, idList.length() - 2);
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 17, idList.c_str());
    }

  if (!this->mFastReactions.empty())
    {
      // create the error message
      std::string idList;
      std::set<std::string>::const_iterator it = this->mFastReactions.begin();
      std::set<std::string>::const_iterator endit = this->mFastReactions.end();

      while (it != endit)
        {
          idList += (*it);
          idList += ", ";
          ++it;
        }

      idList = idList.substr(0, idList.length() - 2);
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 29, idList.c_str());
    }

  if (!this->mReactionsWithReplacedLocalParameters.empty())
    {
      // create the error message
      std::string idList;
      std::set<std::string>::const_iterator it = this->mReactionsWithReplacedLocalParameters.begin();
      std::set<std::string>::const_iterator endit = this->mReactionsWithReplacedLocalParameters.end();

      while (it != endit)
        {
          idList += (*it);
          idList += ", ";
          ++it;
        }

      idList = idList.substr(0, idList.length() - 2);
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 87, idList.c_str());
    }

  // import the initial assignments
  // we do this after the reactions since intial assignments can reference reaction ids.
  if (mpImportHandler)
    {
      mpImportHandler->finishItem(hStep);
      mImportStep = 7;

      if (!mpImportHandler->progressItem(mhImportStep)) return NULL;
    }

  importInitialAssignments(sbmlModel, copasi2sbmlmap, this->mpCopasiModel);

  // import all rules
  // we have to import them after the reactions since they can reference a reaction
  // id which has to be known when importing the rule

  /* Create the rules */
  // rules should be imported after reactions because we use the mSBMLSpeciesReferenceIds to determine if a
  // rule changes a species reference (stoichiometry
  // Since COPASI does not support this, we need to ignore the rule
  this->areRulesUnique(sbmlModel);
  num = sbmlModel->getNumRules();

  if (mpImportHandler)
    {
      mpImportHandler->finishItem(hStep);
      mImportStep = 8;

      if (!mpImportHandler->progressItem(mhImportStep)) return NULL;

      step = 0;
      totalSteps = (unsigned C_INT32) num;
      hStep = mpImportHandler->addItem("Importing global parameters...",
                                       step,
                                       &totalSteps);
    }

  for (counter = 0; counter < num; ++counter)
    {
      Rule* sbmlRule = sbmlModel->getRule(counter);

      if (sbmlRule == NULL)
        {
          fatalError();
        }

      // improve the error message a bit.
      try
        {
          this->importSBMLRule(sbmlRule, copasi2sbmlmap, sbmlModel);
        }
      catch (...)
        {
          std::ostringstream os;
          os << "Error while importing rule for variable \"";
          os << sbmlRule->getVariable() << "\".";

          // check if the last message on the stack is an exception
          // and if so, add the message text to the current exception
          if (CCopasiMessage::peekLastMessage().getType() == CCopasiMessage::EXCEPTION)
            {
              // we only want the message, not the time stamp line
              std::string text = CCopasiMessage::peekLastMessage().getText();
              os << text.substr(text.find("\n"));
            }

          CCopasiMessage(CCopasiMessage::EXCEPTION, os.str().c_str());
        }

      ++step;

      if (mpImportHandler && !mpImportHandler->progressItem(hStep)) return NULL;
    }

  if (sbmlModel->getNumConstraints() > 0)
    {
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 49);
    }

  // TODO Create all constraints
  // TODO Since we don't have constraints yet, there is no code here.
  // TODO When implementing import of constraints, don't forget to replace calls to
  // TODO explicitly time dependent functions in the constraints math expression.

  // import all events
  // events should be imported after reactions because we use the mSBMLSpeciesReferenceIds to determine if an
  // event assignment changes a species reference (stoichiometry
  // Since COPASI does not support this, we need to ignore the event assignment
  if (mpImportHandler)
    {
      mpImportHandler->finishItem(hStep);
      mImportStep = 9;

      if (!mpImportHandler->progressItem(mhImportStep)) return NULL;
    }

  this->importEvents(sbmlModel, this->mpCopasiModel, copasi2sbmlmap);

  this->mpCopasiModel->setCompileFlag();

  if (this->mUnsupportedRuleFound)
    {
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 3);
    }

  if (this->mRateRuleForSpeciesReferenceIgnored == true)
    {
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 94 , "Rate rule" , "Rate rule");
    }

  if (this->mEventAssignmentForSpeciesReferenceIgnored == true)
    {
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 94 , "Event assignment", "event assignment");
    }

  if (this->mConversionFactorFound == true)
    {
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 100);
    }

  if (this->mEventPrioritiesIgnored == true)
    {
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 98);
    }

  if (this->mInitialTriggerValues == true)
    {
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 97);
    }

  if (this->mNonPersistentTriggerFound == true)
    {
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 99);
    }

  if (mpImportHandler)
    {
      mpImportHandler->finishItem(hStep);
      mImportStep = 10;

      if (!mpImportHandler->progressItem(mhImportStep)) return NULL;
    }

  // unset the hasOnlySubstanceUnits flag on all such species
  std::map<Species*, Compartment*>::iterator it = this->mSubstanceOnlySpecies.begin();
  std::map<Species*, Compartment*>::iterator endIt = this->mSubstanceOnlySpecies.end();

  while (it != endIt)
    {
      it->first->setHasOnlySubstanceUnits(false);
      ++it;
    }

  setInitialValues(this->mpCopasiModel, copasi2sbmlmap);
  // evaluate and apply the initial expressions
  this->applyStoichiometricExpressions(copasi2sbmlmap, sbmlModel);

  // now we apply the conversion factors
  if (this->mLevel > 2)
    {
      this->applyConversionFactors();
    }

  this->removeUnusedFunctions(pTmpFunctionDB, copasi2sbmlmap);

  // remove the temporary avogadro parameter if one was created
  if (this->mAvogadroCreated == true)
    {
      const Parameter* pParameter = *this->mPotentialAvogadroNumbers.begin();
      ListOf* pList = sbmlModel->getListOfParameters();
      unsigned i, iMax = pList->size();

      for (i = 0; i < iMax; ++i)
        {
          if (pList->get(i)->getId() == pParameter->getId())
            {
              pList->remove(i);
              break;
            }
        }
    }

  delete pTmpFunctionDB;

  // create a warning if the delay function is used in the model
  if (this->mDelayFound)
    {
      CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 36);
    }

  // give a warning that units on pure number as allowed in SBML L3 and above
  // are ignored by COPASI
  if (this->mLevel > 2 && this->mUnitOnNumberFound)
    {
      CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 93);
    }

  this->mpCopasiModel->forceCompile(this->mpImportHandler);
  return this->mpCopasiModel;
}

CModelValue * SBMLImporter::createCModelValueFromParameter ( const Parameter *  sbmlParameter, CModel *  copasiModel, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap  )

protected

Creates and returns a COPASI CModelValue from the given SBML Parameter object.

Definition at line 4436 of file SBMLImporter.cpp.

References C_INVALID_INDEX, CModel::createModelValue(), CCopasiVectorN< CType >::getIndex(), CModel::getModelValues(), importMIRIAM(), importNotes(), mIgnoredParameterUnits, mLevel, and CModelEntity::setSBMLId().

Referenced by createCModelFromSBMLDocument(), createHasOnlySubstanceUnitFactor(), find_local_parameters_in_delay(), and replace_delay_nodes().

{
  if (sbmlParameter->isSetUnits())
    {
      /* !!!! create a warning that the units will be ignored. */
      //CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 26, sbmlParameter->getId().c_str());
      mIgnoredParameterUnits.push_back(sbmlParameter->getId());
      const_cast<Parameter*>(sbmlParameter)->unsetUnits();
    }

  std::string name = sbmlParameter->getName();

  if (!sbmlParameter->isSetName())
    {
      name = sbmlParameter->getId();
    }

  std::string appendix = "";
  unsigned int counter = 2;
  std::ostringstream numberStream;

  while (copasiModel->getModelValues().getIndex(name + appendix) != C_INVALID_INDEX)
    {
      numberStream.str("");
      numberStream << "_" << counter;
      counter++;
      appendix = numberStream.str();
    }

  std::string sbmlId;

  if (this->mLevel == 1)
    {
      sbmlId = sbmlParameter->getName();
    }
  else
    {
      sbmlId = sbmlParameter->getId();
    }

  CModelValue* pMV = copasiModel->createModelValue(name + appendix, 0.0);
  copasi2sbmlmap[pMV] = const_cast<Parameter*>(sbmlParameter);
  pMV->setSBMLId(sbmlId);
  SBMLImporter::importMIRIAM(sbmlParameter, pMV);
  SBMLImporter::importNotes(pMV, sbmlParameter);

  if (this->mLevel > 2)
    {
      this->mSBMLIdModelValueMap[sbmlId] = pMV;
    }

  return pMV;
}

CReaction * SBMLImporter::createCReactionFromReaction ( Reaction *  sbmlReaction, Model *  pSBMLModel, CModel *  copasiModel, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap, CFunctionDB *  pTmpFunctionDB  )

protected

Creates and returns a COPASI CReaction object from the given SBML Reaction object.

Definition at line 2146 of file SBMLImporter.cpp.

References CFunctionDB::add(), CEvaluationNodeCall::addChild(), CReaction::addModifier(), CCopasiParameterGroup::addParameter(), CReaction::addProduct(), CReaction::addSubstrate(), CCopasiVector< T >::begin(), C_FLOAT64, C_INVALID_INDEX, containsVolume(), CEvaluationNode::copyBranch(), CModel::createReaction(), doMapping(), CCopasiParameter::DOUBLE, CCopasiVector< T >::end(), CCopasiMessage::EXCEPTION, CEvaluationNodeCall::EXPRESSION, fatalError, findCorrespondingFunction(), CFunctionDB::findFunction(), findIdInASTTree(), CEvaluationTree::fromAST(), functionDB, CReaction::getChemEq(), CCopasiNode< _Data >::getChild(), CCopasiObject::getCN(), CCopasiNode< _Data >::getData(), CCompartment::getDimensionality(), CCopasiVectorN< CType >::getIndex(), CChemEq::getModifiers(), CCopasiObject::getObjectName(), CCopasiObject::getObjectType(), CReaction::getParameters(), CChemEq::getProducts(), CModel::getReactions(), CEvaluationTree::getRoot(), CModelEntity::getSBMLId(), CCopasiNode< _Data >::getSibling(), CChemEq::getSubstrates(), CEvaluationTree::getType(), CEvaluationNode::getType(), CFunction::getVariables(), importMIRIAM(), importNotes(), isConstantFlux(), isMassAction(), isMultipliedByVolume(), CReaction::isReversible(), isSimpleFunctionCall(), CFunction::isSuitable(), MCSBML, mDivisionByCompartmentReactions, mFastReactions, mIncompleteModel, mLevel, mStoichiometricExpressionMap, mUsedFunctions, CEvaluationNode::OBJECT, pdelete, preprocessNode(), CChemEq::PRODUCT, renameMassActionParameters(), sbmlId2CopasiCN(), CReaction::setFunction(), CReaction::setFunctionFromExpressionTree(), CCopasiObject::setObjectName(), CFunction::setReversible(), CReaction::setReversible(), CEvaluationTree::setRoot(), CReaction::setSBMLId(), CFunctionParameters::size(), CCopasiVector< T >::size(), speciesMap, CChemEq::SUBSTRATE, TriFalse, TriTrue, TriUnspecified, CEvaluationNode::type(), CEvaluationTree::UserDefined, variables2objects(), and CCopasiMessage::WARNING.

Referenced by createCModelFromSBMLDocument().

{
  if (sbmlReaction == NULL)
    {
      fatalError();
    }

  std::string name = sbmlReaction->getName();

  if (name == "")
    {
      name = sbmlReaction->getId();
    }

  /* Check if the name of the reaction is unique. */
  std::string appendix = "";
  unsigned int counter = 2;
  std::ostringstream numberStream;

  while (copasiModel->getReactions().getIndex(name + appendix) != C_INVALID_INDEX)
    {
      numberStream.str("");
      numberStream << "_" << counter;
      counter++;
      appendix = numberStream.str();
    }

  /* create a new reaction with the unique name */
  CReaction* copasiReaction = copasiModel->createReaction(name + appendix);

  if (copasiReaction == NULL)
    {
      fatalError();
    }

  copasiReaction->setReversible(sbmlReaction->getReversible());

  if (this->mLevel == 1)
    {
      copasiReaction->setSBMLId(sbmlReaction->getName());
    }
  else
    {
      copasiReaction->setSBMLId(sbmlReaction->getId());
    }

  copasi2sbmlmap[copasiReaction] = const_cast<Reaction*>(sbmlReaction);

  if (sbmlReaction->isSetFast() && sbmlReaction->getFast() == true)
    {
      const_cast<Reaction*>(sbmlReaction)->setFast(false);
      this->mFastReactions.insert(sbmlReaction->getId());
    }

  // store if the reaction involves species that need a conversion factor
  bool mConversionFactorNeeded = false;

  /* Add all substrates to the reaction */
  unsigned int num = sbmlReaction->getNumReactants();
  bool singleCompartment = true;
  const CCompartment* compartment = NULL;
  bool hasOnlySubstanceUnitPresent = false;
  bool ignoreMassAction = false;

  for (counter = 0; counter < num; counter++)
    {
      const SpeciesReference* sr = sbmlReaction->getReactant(counter);

      if (sr == NULL)
        {
          delete copasiReaction;
          fatalError();
        }

      C_FLOAT64 stoi = 1.0;

      if (this->mLevel < 3)
        {
          // We may not use Mass Action kinetics if we do not know the actual stoichiometry.
          if (sr->isSetStoichiometryMath())
            {
              ignoreMassAction = true;
            }
          else
            {
              stoi = sr->getStoichiometry() / sr->getDenominator();
            }
        }
      else
        {
          // We may not use Mass Action kinetics if we do not know the actual stoichiometry.
          if (sr->isSetId() &&
              (pSBMLModel->getInitialAssignment(sr->getId()) ||
               pSBMLModel->getRule(sr->getId())))
            {
              ignoreMassAction = true;
            }

          if (sr->isSetStoichiometry())
            {
              stoi = sr->getStoichiometry();
            }
        }

      std::map<std::string, CMetab*>::iterator pos;
      pos = this->speciesMap.find(sr->getSpecies());

      // check if we have a conversion factor
      if (pos == this->speciesMap.end())
        {
          delete copasiReaction;
          fatalError();
        }
      else
        {
          // check if there is a conversion factor on the species
          if (this->mpModelConversionFactor != NULL || this->mSpeciesConversionParameterMap.find(pos->second->getSBMLId()) != this->mSpeciesConversionParameterMap.end())
            {
              mConversionFactorNeeded = true;
            }
        }

      std::map<CCopasiObject*, SBase*>::const_iterator spos = copasi2sbmlmap.find(pos->second);
      assert(spos != copasi2sbmlmap.end());
      Species* pSBMLSpecies = dynamic_cast<Species*>(spos->second);
      assert(pSBMLSpecies != NULL);
      hasOnlySubstanceUnitPresent = (hasOnlySubstanceUnitPresent | (pSBMLSpecies->getHasOnlySubstanceUnits() == true));

      if (compartment == NULL)
        {
          compartment = pos->second->getCompartment();
        }
      else
        {
          if (singleCompartment && compartment != pos->second->getCompartment())
            {
              singleCompartment = false;
            }
        }

      copasiReaction->addSubstrate(pos->second->getKey(), stoi);

      // we need to store the id of the species reference if it is set because SBML Level 3 allows
      // references to species references and if we want to support his, we need the id to import
      // expressions that reference a species reference
      if (this->mLevel > 2)
        {
          CCopasiVector<CChemEqElement>::const_iterator it = copasiReaction->getChemEq().getSubstrates().begin(), endit = copasiReaction->getChemEq().getSubstrates().end();
          CChemEqElement* pElement = NULL;

          while (it != endit)
            {
              if ((*it)->getMetaboliteKey() == pos->second->getKey())
                {
                  pElement = const_cast<CChemEqElement*>(*it);
                  break;
                }

              ++it;
            }

          assert(pElement != NULL);

          if (pElement != NULL)
            {
              copasi2sbmlmap[pElement] = const_cast<SpeciesReference*>(sr);

              this->mChemEqElementSpeciesIdMap[pElement] = std::pair<std::string, CChemEq::MetaboliteRole>(sr->getSpecies(), CChemEq::SUBSTRATE);
            }
        }

      // find the CChemEqElement that belongs to the added substrate
      if (this->mLevel < 3 && sr->isSetStoichiometryMath())
        {
          CChemEq& chemEq = copasiReaction->getChemEq();
          CCopasiVector < CChemEqElement >::const_iterator it = chemEq.getSubstrates().begin();
          CCopasiVector < CChemEqElement >::const_iterator end = chemEq.getSubstrates().end();
          CChemEqElement* pChemEqElement = NULL;

          while (it != end)
            {
              if ((*it)->getMetabolite() == pos->second)
                {
                  pChemEqElement = (*it);
                  break;
                }

              ++it;
            }

          assert(pChemEqElement != NULL);
          mStoichiometricExpressionMap.insert(std::make_pair(sr->getStoichiometryMath()->getMath(), pChemEqElement));
        }
    }

  /* Add all products to the reaction */
  num = sbmlReaction->getNumProducts();

  for (counter = 0; counter < num; counter++)
    {
      const SpeciesReference* sr = sbmlReaction->getProduct(counter);

      if (sr == NULL)
        {
          delete copasiReaction;
          fatalError();
        }

      C_FLOAT64 stoi = 1.0;

      if (this->mLevel < 3)
        {
          if (sr->isSetStoichiometryMath())
            {
              ignoreMassAction = true;
            }
          else
            {
              stoi = sr->getStoichiometry() / sr->getDenominator();
            }
        }
      else
        {
          // We may not use Mass Action kinetics if we do not know the actual stoichiometry.
          if (sr->isSetId() &&
              (pSBMLModel->getInitialAssignment(sr->getId()) ||
               pSBMLModel->getRule(sr->getId())))
            {
              ignoreMassAction = true;
            }

          if (sr->isSetStoichiometry())
            {
              stoi = sr->getStoichiometry();
            }
        }

      std::map<std::string, CMetab*>::iterator pos;
      pos = this->speciesMap.find(sr->getSpecies());

      if (pos == this->speciesMap.end())
        {
          delete copasiReaction;
          fatalError();
        }
      else
        {
          // check if there is a conversion factor on the species
          if (this->mpModelConversionFactor != NULL || this->mSpeciesConversionParameterMap.find(pos->second->getSBMLId()) != this->mSpeciesConversionParameterMap.end())
            {
              mConversionFactorNeeded = true;
            }
        }

      std::map<CCopasiObject*, SBase*>::const_iterator spos = copasi2sbmlmap.find(pos->second);
      assert(spos != copasi2sbmlmap.end());
      Species* pSBMLSpecies = dynamic_cast<Species*>(spos->second);
      assert(pSBMLSpecies != NULL);
      hasOnlySubstanceUnitPresent = (hasOnlySubstanceUnitPresent | (pSBMLSpecies->getHasOnlySubstanceUnits() == true));

      if (compartment == NULL)
        {
          compartment = pos->second->getCompartment();
        }
      else
        {
          if (singleCompartment && compartment != pos->second->getCompartment())
            {
              singleCompartment = false;
            }
        }

      copasiReaction->addProduct(pos->second->getKey(), stoi);

      // we need to store the id of the species reference if it is set because SBML Level 3 allows
      // references to species references and if we want to support his, we need the id to import
      // expressions that reference a species reference
      if (this->mLevel > 2)
        {
          CCopasiVector<CChemEqElement>::const_iterator it = copasiReaction->getChemEq().getProducts().begin(), endit = copasiReaction->getChemEq().getProducts().end();
          CChemEqElement* pElement = NULL;

          while (it != endit)
            {
              if ((*it)->getMetaboliteKey() == pos->second->getKey())
                {
                  pElement = const_cast<CChemEqElement*>(*it);
                  break;
                }

              ++it;
            }

          assert(pElement != NULL);

          if (pElement != NULL)
            {
              copasi2sbmlmap[pElement] = const_cast<SpeciesReference*>(sr);
              this->mChemEqElementSpeciesIdMap[pElement] = std::pair<std::string, CChemEq::MetaboliteRole>(sr->getSpecies(), CChemEq::PRODUCT);
            }
        }

      if (sr->isSetStoichiometryMath())
        {
          CChemEq& chemEq = copasiReaction->getChemEq();
          CCopasiVector < CChemEqElement >::const_iterator it = chemEq.getProducts().begin();
          CCopasiVector < CChemEqElement >::const_iterator end = chemEq.getProducts().end();
          CChemEqElement* pChemEqElement = NULL;

          while (it != end)
            {
              if ((*it)->getMetabolite() == pos->second)
                {
                  pChemEqElement = (*it);
                  break;
                }

              ++it;
            }

          assert(pChemEqElement != NULL);
          mStoichiometricExpressionMap.insert(std::make_pair(sr->getStoichiometryMath()->getMath(), pChemEqElement));
        }
    }

  /* Add all modifiers to the reaction */
  num = sbmlReaction->getNumModifiers();

  for (counter = 0; counter < num; counter++)
    {
      const ModifierSpeciesReference* sr = sbmlReaction->getModifier(counter);

      if (sr == NULL)
        {
          delete copasiReaction;
          fatalError();
        }

      std::map<std::string, CMetab*>::iterator pos;
      pos = this->speciesMap.find(sr->getSpecies());

      if (pos == this->speciesMap.end())
        {
          delete copasiReaction;
          fatalError();
        }

      std::map<CCopasiObject*, SBase*>::const_iterator spos = copasi2sbmlmap.find(pos->second);
      assert(spos != copasi2sbmlmap.end());
      Species* pSBMLSpecies = dynamic_cast<Species*>(spos->second);
      assert(pSBMLSpecies != NULL);
      hasOnlySubstanceUnitPresent = (hasOnlySubstanceUnitPresent | (pSBMLSpecies->getHasOnlySubstanceUnits() == true));
      copasiReaction->addModifier(pos->second->getKey());

      // we need to store the id of the species reference if it is set because SBML Level 3 allows
      // references to species references and if we want to support his, we need the id to import
      // expressions that reference a species reference
      if (this->mLevel > 2)
        {

          CCopasiVector<CChemEqElement>::const_iterator it = copasiReaction->getChemEq().getModifiers().begin(), endit = copasiReaction->getChemEq().getModifiers().end();
          CChemEqElement* pElement = NULL;

          while (it != endit)
            {
              if ((*it)->getMetaboliteKey() == pos->second->getKey())
                {
                  pElement = const_cast<CChemEqElement*>(*it);
                  break;
                }

              ++it;
            }

          assert(pElement != NULL);

          if (pElement != NULL)
            {
              copasi2sbmlmap[pElement] = const_cast<ModifierSpeciesReference*>(sr);
            }
        }
    }

  /* in the newly created CFunction set the types for all parameters and
   * either a mapping or a value
   */
  const KineticLaw* kLaw = sbmlReaction->getKineticLaw();

  if (kLaw != NULL)
    {
      const ListOfParameters* pParamList = NULL;

      if (this->mLevel > 2)
        {
          pParamList = kLaw->getListOfLocalParameters();
        }
      else
        {
          pParamList = kLaw->getListOfParameters();
        }

      for (counter = 0; counter < pParamList->size(); ++counter)
        {
          const Parameter* pSBMLParameter = pParamList->get(counter);
          std::string id;

          if (this->mLevel == 1)
            {
              id = pSBMLParameter->getName();
            }
          else
            {
              id = pSBMLParameter->getId();
            }

          double value;

          if (pSBMLParameter->isSetValue() && pSBMLParameter->getValue() == pSBMLParameter->getValue()) // make sure it is not set to NaN
            {
              value = pSBMLParameter->getValue();
            }
          else
            {
              // Set value to NaN and create a warning if it is the first time
              // this happend
              value = std::numeric_limits<C_FLOAT64>::quiet_NaN();

              if (!this->mIncompleteModel)
                {
                  this->mIncompleteModel = true;
                  CCopasiMessage Message(CCopasiMessage::WARNING, MCSBML + 42, pSBMLParameter->getId().c_str());
                }
            }

          copasiReaction->getParameters().addParameter(id, CCopasiParameter::DOUBLE, value);
        }

      const ASTNode* kLawMath = kLaw->getMath();

      if (kLawMath == NULL || kLawMath->getType() == AST_UNKNOWN)
        {
          copasiReaction->setFunction(NULL);
          CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 56, sbmlReaction->getId().c_str());
        }
      else
        {
          // check for references to species references in the expression because we don't support them yet
          if (!SBMLImporter::findIdInASTTree(kLawMath, this->mSBMLSpeciesReferenceIds).empty())
            {
              CCopasiMessage(CCopasiMessage::WARNING, MCSBML + 95);
            }

          ConverterASTNode* node = new ConverterASTNode(*kLawMath);
          this->preprocessNode(node, pSBMLModel, copasi2sbmlmap, sbmlReaction);

          if (node == NULL)
            {
              delete copasiReaction;
              fatalError();
            }

          /* if it is a single compartment reaction, we have to divide the whole kinetic
          ** equation by the compartment because COPASI expects
          ** kinetic laws that specify concentration/time for single compartment
          ** reactions.
          */
          if (singleCompartment)
            {

              if (compartment != NULL)
                {
                  // only divide if it is not a 0-dimensional compartment
                  if (compartment->getDimensionality() != 0)
                    {
                      // check if the root node is a multiplication node and it's first child
                      // is a compartment node, if so, drop those two and make the second child
                      // new new root node
                      ConverterASTNode* tmpNode1 = this->isMultipliedByVolume(node, compartment->getSBMLId());

                      if (tmpNode1)
                        {
                          delete node;
                          node = tmpNode1;

                          if (node->getType() == AST_DIVIDE && node->getNumChildren() != 2)
                            {
                              delete tmpNode1;
                              fatalError();
                            }
                        }
                      else
                        {
                          tmpNode1 = new ConverterASTNode();
                          tmpNode1->setType(AST_DIVIDE);
                          tmpNode1->addChild(node);
                          ConverterASTNode* tmpNode2 = new ConverterASTNode();
                          tmpNode2->setType(AST_NAME);
                          tmpNode2->setName(compartment->getSBMLId().c_str());
                          tmpNode1->addChild(tmpNode2);
                          node = tmpNode1;
                          std::map<CCopasiObject*, SBase*>::const_iterator pos = copasi2sbmlmap.find(const_cast<CCompartment*>(compartment));
                          assert(pos != copasi2sbmlmap.end());
                          Compartment* pSBMLCompartment = dynamic_cast<Compartment*>(pos->second);
                          assert(pSBMLCompartment != NULL);

                          if (!hasOnlySubstanceUnitPresent && ((this->mLevel == 1 && pSBMLCompartment->isSetVolume()) || (this->mLevel >= 2 && pSBMLCompartment->isSetSize())) && pSBMLCompartment->getSize() == 1.0)
                            {
                              // we have to check if all species used in the reaction
                              // have the hasOnlySubstance flag set

                              if (node->getChild(0)->getType() == AST_FUNCTION && (!this->containsVolume(node->getChild(0), compartment->getSBMLId())))
                                {
                                  // add the id of the reaction to the set so that we can create an error message later.
                                  this->mDivisionByCompartmentReactions.insert(sbmlReaction->getId());
                                }
                            }
                        }
                    }
                }
              else
                {
                  delete node;
                  delete copasiReaction;
                  fatalError();
                }
            }

          /* Create a new user defined CKinFunction */
          if (!sbmlId2CopasiCN(node, copasi2sbmlmap, copasiReaction->getParameters()))
            {
              CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 27, copasiReaction->getObjectName().c_str());
            }

          CEvaluationNode* pExpressionTreeRoot = CEvaluationTree::fromAST(node);
          delete node;
          node = NULL;

          if (pExpressionTreeRoot)
            {
              CExpression KineticLawExpression;
              KineticLawExpression.setRoot(pExpressionTreeRoot);

              // check if the expression is constant flux
              CCopasiObject* pParamObject = SBMLImporter::isConstantFlux(pExpressionTreeRoot, copasiModel, pTmpFunctionDB);

              if (pParamObject != NULL)
                {
                  // assert that the object really is a local or global
                  // parameter
                  assert(dynamic_cast<const CCopasiParameter*>(pParamObject) || dynamic_cast<const CModelValue*>(pParamObject));
                  std::string functionName;

                  if (copasiReaction->isReversible())
                    {
                      // set the function to Constant flux (reversible)
                      functionName = "Constant flux (reversible)";
                    }
                  else
                    {
                      // set the function to Constant flux (irreversible)
                      functionName = "Constant flux (irreversible)";
                    }

                  CFunction* pCFFun = dynamic_cast<CFunction*>(this->functionDB->findFunction(functionName));
                  assert(pCFFun != NULL);
                  CEvaluationNodeCall* pCallNode = NULL;

                  if (CEvaluationNode::type(pExpressionTreeRoot->getType()) == CEvaluationNode::OBJECT)
                    {
                      pCallNode = new CEvaluationNodeCall(CEvaluationNodeCall::EXPRESSION, "dummy_call");
                      // add the parameter
                      pCallNode->addChild(pExpressionTreeRoot->copyBranch());
                    }
                  else
                    {
                      pCallNode = dynamic_cast<CEvaluationNodeCall*>(pExpressionTreeRoot->copyBranch());
                      assert(pCallNode != NULL);
                    }

                  if (pParamObject->getObjectType() == "Parameter")
                    {
                      pParamObject->setObjectName("v");
                      dynamic_cast<CEvaluationNode*>(pCallNode->getChild())->setData("<" + pParamObject->getCN() + ">");
                    }

                  copasiReaction->setFunction(pCFFun);
                  // map the parameter
                  this->doMapping(copasiReaction, pCallNode);
                  delete pCallNode;
                }
              else
                {
                  // check if the root node is a simple function call
                  if (this->isSimpleFunctionCall(pExpressionTreeRoot))
                    {
                      // if yes, we check if it corresponds to an already existing function
                      std::string functionName = pExpressionTreeRoot->getData();
                      CFunction* pImportedFunction = dynamic_cast<CFunction*>(functionDB->findFunction(functionName));
                      assert(pImportedFunction);
                      std::vector<CEvaluationNodeObject*>* v = NULL;

                      // only check for mass action if there is no conversion factor involved
                      // for any of the species involved in the reaction (substrates and products)
                      if (!mConversionFactorNeeded && !ignoreMassAction)
                        {
                          v = this->isMassAction(pImportedFunction, copasiReaction->getChemEq(), static_cast<const CEvaluationNodeCall*>(pExpressionTreeRoot));

                          if (!v)
                            {
                              CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 27, copasiReaction->getObjectName().c_str());
                            }
                        }

                      if (v && !v->empty())
                        {
                          CFunction* pFun = NULL;

                          if (copasiReaction->isReversible())
                            {
                              pFun = static_cast<CFunction*>(this->functionDB->findFunction("Mass action (reversible)"));
                            }
                          else
                            {
                              pFun = static_cast<CFunction*>(this->functionDB->findFunction("Mass action (irreversible)"));
                            }

                          if (!pFun)
                            {
                              fatalError();
                            }

                          // do the mapping
                          CEvaluationNodeCall* pCallNode = new CEvaluationNodeCall(CEvaluationNodeCall::EXPRESSION, "dummy_call");
                          size_t i, iMax = v->size();

                          for (i = 0; i < iMax; ++i)
                            {
                              pCallNode->addChild((*v)[i]);
                            }

                          this->renameMassActionParameters(pCallNode);
                          copasiReaction->setFunction(pFun);
                          this->doMapping(copasiReaction, pCallNode);
                          delete pCallNode;
                        }
                      else
                        {
                          // find corresponding reaction does *not* return the function if it is already in the function db
                          // we better change this to return this instance.
                          CFunction* pExistingFunction = findCorrespondingFunction(&KineticLawExpression, copasiReaction);

                          // if it does, we set the existing function for this reaction
                          if (pExistingFunction)
                            {
                              copasiReaction->setFunction(pExistingFunction);
                              // do the mapping
                              doMapping(copasiReaction, dynamic_cast<const CEvaluationNodeCall*>(pExpressionTreeRoot));
                            }
                          // else we take the function from the pTmpFunctionDB, copy it and set the usage correctly
                          else
                            {
                              // replace the variable nodes in pImportedFunction with  nodes from
                              std::map<std::string, std::string > arguments;
                              const CFunctionParameters& funParams = pImportedFunction->getVariables();
                              const CEvaluationNode* pTmpNode = static_cast<const CEvaluationNode*>(pExpressionTreeRoot->getChild());
                              size_t i, iMax = funParams.size();

                              for (i = 0; (i < iMax) && pTmpNode; ++i)
                                {
                                  if (!(pTmpNode->getType() == CEvaluationNode::OBJECT)) fatalError();

                                  arguments[funParams[i]->getObjectName()] = pTmpNode->getData().substr(1, pTmpNode->getData().length() - 2);
                                  pTmpNode = static_cast<const CEvaluationNode*>(pTmpNode->getSibling());
                                }

                              assert((i == iMax) && pTmpNode == NULL);

                              CEvaluationNode* pTmpExpression = variables2objects(pImportedFunction->getRoot()->copyBranch(), arguments);

                              CExpression TmpTree2;
                              TmpTree2.setRoot(pTmpExpression);

                              // code to fix bug 1874
                              // since this is a user defined function, we want to retain the original name
                              // next setFunctionFromExpressionTree will take this name if it is not 'Expression'
                              // (the default)
                              TmpTree2.setObjectName(functionName);

                              CFunction * pNewFunction = copasiReaction->setFunctionFromExpressionTree(TmpTree2, copasi2sbmlmap, this->functionDB);

                              if (pNewFunction != NULL &&
                                  pNewFunction->getType() == CEvaluationTree::UserDefined)
                                {
                                  // code to fix Bug 1015
                                  if (!pNewFunction->isSuitable(copasiReaction->getChemEq().getSubstrates().size(),
                                                                copasiReaction->getChemEq().getProducts().size(),
                                                                copasiReaction->isReversible() ? TriTrue : TriFalse))
                                    {
                                      pNewFunction->setReversible(TriUnspecified);
                                    }

                                  pTmpFunctionDB->add(pNewFunction, false);
                                  mUsedFunctions.insert(pNewFunction->getObjectName());
                                }
                            }
                        }

                      pdelete(v);
                    }
                  else
                    {
                      std::vector<CEvaluationNodeObject*>* v = NULL;

                      // only check for mass action if there is no conversion factor involved
                      // for any of the species involved in the reaction (substrates and products)
                      if (!mConversionFactorNeeded && !ignoreMassAction)
                        {
                          v = this->isMassAction(&KineticLawExpression, copasiReaction->getChemEq());

                          if (!v)
                            {
                              CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 27, copasiReaction->getObjectName().c_str());
                            }
                        }

                      if (v && !v->empty())
                        {
                          CFunction* pFun = NULL;

                          if (copasiReaction->isReversible())
                            {
                              pFun = static_cast<CFunction*>(this->functionDB->findFunction("Mass action (reversible)"));
                            }
                          else
                            {
                              pFun = static_cast<CFunction*>(this->functionDB->findFunction("Mass action (irreversible)"));
                            }

                          if (!pFun)
                            {
                              fatalError();
                            }

                          // do the mapping
                          CEvaluationNodeCall* pCallNode = new CEvaluationNodeCall(CEvaluationNodeCall::EXPRESSION, "dummy_call");
                          size_t i, iMax = v->size();

                          for (i = 0; i < iMax; ++i)
                            {
                              pCallNode->addChild((*v)[i]);
                            }

                          // rename the function parameters to k1 and k2
                          this->renameMassActionParameters(pCallNode);
                          copasiReaction->setFunction(pFun);
                          this->doMapping(copasiReaction, pCallNode);
                          delete pCallNode;
                        }
                      else
                        {
                          CFunction* pNonconstFun = copasiReaction->setFunctionFromExpressionTree(KineticLawExpression, copasi2sbmlmap, this->functionDB);

                          if (pNonconstFun != NULL &&
                              pNonconstFun->getType() == CEvaluationTree::UserDefined)
                            {
                              // code to fix Bug 1015
                              if (!pNonconstFun->isSuitable(copasiReaction->getChemEq().getSubstrates().size(),
                                                            copasiReaction->getChemEq().getProducts().size(),
                                                            copasiReaction->isReversible() ? TriTrue : TriFalse))
                                {
                                  pNonconstFun->setReversible(TriUnspecified);
                                }

                              pTmpFunctionDB->add(pNonconstFun, false);
                              mUsedFunctions.insert(pNonconstFun->getObjectName());
                            }
                        }

                      pdelete(v);
                    }
                }
            }
          else
            {
              // error message
              CCopasiMessage(CCopasiMessage::EXCEPTION, MCSBML + 8, copasiReaction->getObjectName().c_str());
            }
        }
    }
  else
    {
      /* if no KineticLaw was defined for the reaction. */
      copasiReaction->setFunction(NULL);
    }

  //DebugFile << "Created reaction: " << copasiReaction->getObjectName() << std::endl;
  SBMLImporter::importMIRIAM(sbmlReaction, copasiReaction);
  SBMLImporter::importNotes(copasiReaction, sbmlReaction);

  return copasiReaction;
}

void SBMLImporter::createDelayFunctionDefinition ( )

protected

Creates a function definition for the delay function.

CEvaluationTree * SBMLImporter::createExpressionFromFunction ( const CFunction *  pFun, const std::vector< std::vector< std::string > > &  functionArgumentCNs  )

protected

Definition at line 5516 of file SBMLImporter.cpp.

References CEvaluationTree::create(), CEvaluationTree::Expression, CCopasiObject::getObjectName(), CEvaluationTree::getRoot(), CFunction::getVariables(), CEvaluationTree::setRoot(), CFunctionParameters::size(), and variables2objects().

Referenced by isMassActionFunction().

{
  CEvaluationTree* pTree = NULL;
  const CFunctionParameters& pFunParams = pFun->getVariables();
  std::string str;

  if (pFunParams.size() == functionArgumentCNs.size())
    {
      std::map<std::string , std::string> variable2CNMap;
      size_t i, iMax = pFunParams.size();

      for (i = 0; i < iMax; ++i)
        {
          // vectors should not occur here
          assert(functionArgumentCNs[i].size() == 1);
          variable2CNMap[pFunParams[i]->getObjectName()] = functionArgumentCNs[i][0];
        }

      CEvaluationNode* pTmpNode = this->variables2objects(pFun->getRoot(), variable2CNMap);
      assert(pTmpNode);
      pTree = CEvaluationTree::create(CEvaluationTree::Expression);
      pTree->setRoot(pTmpNode);
    }

  return pTree;
}

void SBMLImporter::createHasOnlySubstanceUnitFactor ( Model *  pSBMLModel, double  factor, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap  )

protected

This method creates a global parameter the represents the factor that is used to convert a particle number into the amount units set on the model. The parameter is only created if it is needed and after exporting the model, the parameter is deleted from the COPASI model again.

Definition at line 9106 of file SBMLImporter.cpp.

References createCModelValueFromParameter(), mAvogadroCreated, mpCopasiModel, and mPotentialAvogadroNumbers.

Referenced by replaceAmountReferences().

{
  // find an ID that is unique at least within the list of parameters
  // since we remove this created parameter after import, it does not
  // have to be unique within the whole model
  std::set<std::string> ids;
  unsigned int i, iMax = pSBMLModel->getListOfParameters()->size();

  for (i = 0; i < iMax; ++i)
    {
      ids.insert(pSBMLModel->getListOfParameters()->get(i)->getId());
    }

  std::ostringstream os;
  i = 1;
  os << "parameter_" << i;

  while (ids.find(os.str()) != ids.end())
    {
      ++i;
      os.str("");
      os << "parameter_" << i;
    }

  Parameter *pParameter = pSBMLModel->createParameter();
  pParameter->setId(os.str());
  pParameter->setName("amount to particle factor");
  pParameter->setConstant(true);
  pParameter->setValue(factor);
  this->mAvogadroCreated = true;
  this->mPotentialAvogadroNumbers.insert(pParameter);
  this->createCModelValueFromParameter(pParameter, this->mpCopasiModel, copasi2sbmlmap);
}

void SBMLImporter::deleteCopasiModel

(

)

This call deletes an existing COPASI model. The method can e.g. be called to clean up if an import fails.

Definition at line 10681 of file SBMLImporter.cpp.

References mpCopasiModel.

Referenced by SEDMLImporter::importFirstSBMLModel(), CCopasiDataModel::importSBML(), and CCopasiDataModel::importSBMLFromString().

{
  if (this->mpCopasiModel != NULL)
    {
      delete this->mpCopasiModel;
      this->mpCopasiModel = NULL;
    }
}

CEvaluationNode * SBMLImporter::divideByObject ( const CEvaluationNode *  pOrigNode, const CCopasiObject *  pObject  )

static

This method divides the given expression by the given object and returns a new expression. The caller is responsible for freeing the memory for the new expression.

Definition at line 10533 of file SBMLImporter.cpp.

References CCopasiNode< _Data >::addChild(), CEvaluationNodeObject::CN, CEvaluationNode::copyBranch(), CEvaluationNodeOperator::DIVIDE, CCopasiNode< _Data >::getChild(), CCopasiObject::getCN(), CCopasiNode< _Data >::getSibling(), CEvaluationNode::getType(), CEvaluationNodeOperator::MULTIPLY, CEvaluationNode::OBJECT, CEvaluationNode::OPERATOR, pResult, CEvaluationNode::subType(), and CEvaluationNode::type().

Referenced by importEvent(), importInitialAssignments(), and importRuleForModelEntity().

{
  bool reverse = false;
  CEvaluationNode* pResult = NULL;
  assert(pOrigNode != NULL);
  assert(pObject != NULL);

  if (pOrigNode != NULL && pObject != NULL)
    {
      // first we check if this is thie reverse operation with the object
      // if so, we just drop the reverse operaqtion, otherwise we apply the operation with the
      // object
      if (CEvaluationNode::type(pOrigNode->getType()) == CEvaluationNode::OPERATOR &&
          (CEvaluationNodeOperator::SubType)CEvaluationNode::subType(pOrigNode->getType()) == CEvaluationNodeOperator::MULTIPLY)
        {
          // either child can be the object
          const CEvaluationNode* pChild = dynamic_cast<const CEvaluationNode*>(pOrigNode->getChild());

          if (CEvaluationNode::type(pChild->getType()) == CEvaluationNode::OBJECT && dynamic_cast<const CEvaluationNodeObject*>(pChild)->getData() == std::string("<" + pObject->getCN() + ">"))
            {

              pResult = dynamic_cast<const CEvaluationNode*>(pOrigNode->getChild())->copyBranch();
              reverse = true;
            }

          if (reverse == false)
            {
              pChild = dynamic_cast<const CEvaluationNode*>(pChild->getSibling());

              if (CEvaluationNode::type(pChild->getType()) == CEvaluationNode::OBJECT && dynamic_cast<const CEvaluationNodeObject*>(pChild)->getData() == std::string("<" + pObject->getCN() + ">"))
                {

                  pResult = dynamic_cast<const CEvaluationNode*>(pOrigNode->getChild())->copyBranch();
                  reverse = true;
                }
            }
        }

      if (reverse == false)
        {
          CEvaluationNodeObject* pVolumeNode = new CEvaluationNodeObject(CEvaluationNodeObject::CN, "<" + pObject->getCN() + ">");
          pResult = new CEvaluationNodeOperator(CEvaluationNodeOperator::DIVIDE, "/");
          pResult->addChild(pOrigNode->copyBranch());
          pResult->addChild(pVolumeNode);
        }
    }

  return pResult;
}

void SBMLImporter::doMapping ( CReaction *  pCopasiReaction, const CEvaluationNodeCall *  pCallNode  )

protected

Definition at line 5740 of file SBMLImporter.cpp.

References CReaction::addParameterMapping(), containsKey(), fatalError, CReaction::getChemEq(), CCopasiNode< _Data >::getChild(), CEvaluationNodeObject::getData(), CReaction::getFunction(), CCopasiObject::getKey(), CChemEq::getModifiers(), CCopasiObject::getObjectParent(), CCopasiObject::getObjectType(), CChemEq::getProducts(), CChemEq::getReversibility(), CCopasiNode< _Data >::getSibling(), CChemEq::getSubstrates(), CFunction::getVariables(), CCopasiObject::isReference(), CReaction::isReversible(), CFunctionParameter::MODIFIER, mpCopasiModel, mpDataModel, CCopasiDataModel::ObjectFromName(), CFunctionParameter::PARAMETER, CFunctionParameter::PRODUCT, CReaction::setParameterMapping(), CFunctionParameters::size(), CCopasiVector< T >::size(), CFunctionParameter::SUBSTRATE, CFunctionParameter::TIME, CFunctionParameter::VARIABLE, and CFunctionParameter::VOLUME.

Referenced by createCReactionFromReaction().

{
  // map the first argument of the call node to the first variable of the function of the reaction
  // and so on
  if (!pCallNode)
    {
      fatalError();
    }

  std::vector<CCopasiContainer*> listOfContainers;
  listOfContainers.push_back(this->mpCopasiModel);

  if (dynamic_cast<const CMassAction*>(pCopasiReaction->getFunction()))
    {
      const CEvaluationNodeObject* pChild = dynamic_cast<const CEvaluationNodeObject*>(pCallNode->getChild());
      std::string objectCN = pChild->getData();
      objectCN = objectCN.substr(1, objectCN.length() - 2);
      CCopasiObject* pObject = mpDataModel->ObjectFromName(listOfContainers, objectCN);

      if (!pObject)
        {
          fatalError();
        }

      if (pObject->isReference())
        {
          pObject = pObject->getObjectParent();
        }

      const std::string& objectKey = pObject->getKey();

      pCopasiReaction->setParameterMapping("k1", objectKey);

      const CCopasiVector<CChemEqElement>* metabolites = &pCopasiReaction->getChemEq().getSubstrates();

      size_t i, iMax = metabolites->size();

      size_t j, jMax;

      for (i = 0; i < iMax; ++i)
        for (j = 0, jMax = static_cast<int>(fabs((*metabolites)[i]->getMultiplicity())); j < jMax; j++)
          pCopasiReaction->addParameterMapping("substrate", (*metabolites)[i]->getMetaboliteKey());

      if (pCopasiReaction->isReversible())
        {
          pChild = dynamic_cast<const CEvaluationNodeObject*>(pChild->getSibling());
          std::string objectCN = pChild->getData();
          objectCN = objectCN.substr(1, objectCN.length() - 2);
          CCopasiObject* pObject = mpDataModel->ObjectFromName(listOfContainers, objectCN);

          if (!pObject)
            {
              fatalError();
            }

          if (pObject->isReference())
            {
              pObject = pObject->getObjectParent();
            }

          const std::string& objectKey = pObject->getKey();

          pCopasiReaction->setParameterMapping("k2", objectKey);

          const CCopasiVector<CChemEqElement>* metabolites = &pCopasiReaction->getChemEq().getProducts();

          iMax = metabolites->size();

          for (i = 0; i < iMax; ++i)
            for (j = 0, jMax = static_cast<int>(fabs((*metabolites)[i]->getMultiplicity())); j < jMax; j++)
              pCopasiReaction->addParameterMapping("product", (*metabolites)[i]->getMetaboliteKey());
        }
    }
  else
    {
      const CFunctionParameters & Variables = pCopasiReaction->getFunction()->getVariables();
      size_t i, iMax = Variables.size();
      const CEvaluationNodeObject* pChild = dynamic_cast<const CEvaluationNodeObject*>(pCallNode->getChild());

      for (i = 0; i < iMax; ++i)
        {
          if (!pChild)
            {
              fatalError();
            }

          std::string objectCN = pChild->getData();
          objectCN = objectCN.substr(1, objectCN.length() - 2);
          CCopasiObject* pObject = mpDataModel->ObjectFromName(listOfContainers, objectCN);

          if (!pObject)
            {
              fatalError();
            }

          if (pObject->isReference())
            {
              pObject = pObject->getObjectParent();
            }

          const std::string& objectKey = pObject->getKey();

          pCopasiReaction->setParameterMapping(i, objectKey);

          const CChemEq& eqn = pCopasiReaction->getChemEq();

          bool reversible = eqn.getReversibility();

          // We guess what the role of a variable of newly imported function is:
          if (Variables[i]->getUsage() == CFunctionParameter::VARIABLE)
            {
              CFunctionParameter::Role Role = CFunctionParameter::PARAMETER;

              if (pObject->getObjectType() == "Metabolite")
                {
                  if (containsKey(eqn.getSubstrates(), objectKey))
                    Role = CFunctionParameter::SUBSTRATE;

                  if (Role == CFunctionParameter::PARAMETER &&
                      containsKey(eqn.getProducts(), objectKey))
                    {
                      // Fix for Bug 1882: if not reversible, only mark as product if it does
                      // not appear in the list of modifiers
                      if (!reversible && containsKey(eqn.getModifiers(), objectKey))
                        Role =  CFunctionParameter::MODIFIER;
                      else
                        Role =  CFunctionParameter::PRODUCT;
                    }

                  // It is not a substrate and not a product therefore we must have a modifier
                  if (Role == CFunctionParameter::PARAMETER)
                    {
                      Role = CFunctionParameter::MODIFIER;
                    }
                }
              else if (pObject->getObjectType() == "Model")
                {
                  Role = CFunctionParameter::TIME;
                }
              else if (pObject->getObjectType() == "Compartment")
                {
                  Role = CFunctionParameter::VOLUME;
                }

              const_cast< CFunctionParameter * >(Variables[i])->setUsage(Role);
            }

          pChild = dynamic_cast<const CEvaluationNodeObject*>(pChild->getSibling());
        }
    }
}

void SBMLImporter::find_local_parameters_in_delay ( ASTNode *  pASTNode, Reaction *  pSBMLReaction, Model *  pModel, std::map< std::string, std::string > &  localReplacementMap, const std::set< std::string > &  localIds, std::map< CCopasiObject *, SBase * > &  copasi2sbmlmap  )

protected

If we replace delay nodes within kineitc laws, we have to make sure that there is no reference to a local parameter within the replaced delay node because that would mean that we end up with a reference to a local parameter in the rule for the delay replacement which is not allowed in SBML. Therefore we have to convert all local parameters which occur within a delay call into global parameters. This method finds all local parameters that have to be converted to global parameters and it already creates the necessary global parameters. The localReplacementMap returns a mapping between the is of the original parameter and the id of the new parameter it will be replaced with. This map is used in a second step to actually replace the nodes in the expression.

If we replace delay nodes within kinetic laws, we have to make sure that there is no reference to a local parameter within the replaced delay node because that would mean that we end up with a reference to a local parameter in the rule for the delay replacement which is not allowed in SBML. Therefore we have to convert all local parameters which occur within a delay call into global parameters.

Definition at line 10234 of file SBMLImporter.cpp.

References createCModelValueFromParameter(), CNodeContextIterator< Node, int >::end(), CCopasiMessage::EXCEPTION, fatalError, mpCopasiModel, mUsedSBMLIds, and CNodeContextIterator< Node, int >::next().

Referenced by replace_delay_nodes().

{
  CNodeIterator< ASTNode > itNode(pASTNode);

  while (itNode.next() != itNode.end())
    {
      if (*itNode == NULL)
        {
          continue;
        }

      if (itNode->getType() == AST_NAME)
        {
          // check it this is the name of a local parameter

          // maybe there already is a replacement for this
          std::map<std::string, std::string>::const_iterator pos = localReplacementMap.find(itNode->getName());