COPASI API  4.16.103
Public Member Functions | Protected Attributes | List of all members
CSEDMLExporter Class Reference

#include <CSEDMLExporter.h>

Public Member Functions

void createDataGenerators (CCopasiDataModel &dataModel, std::string &taskId, CCopasiTask *task=NULL)
 
void createModels (CCopasiDataModel &dataModel, std::string &modelRef)
 
std::string createScanTask (CCopasiDataModel &dataModel, const std::string &modelId)
 
void createSEDMLDocument (CCopasiDataModel &dataModel, std::string modelRef)
 
std::string createSteadyStateTask (CCopasiDataModel &dataModel, const std::string &modelId)
 
void createTasks (CCopasiDataModel &dataModel, std::string &modelRef)
 
std::string createTimeCourseTask (CCopasiDataModel &dataModel, const std::string &modelId)
 
 CSEDMLExporter ()
 
bool exportModelAndTasks (CCopasiDataModel &dataModel, const std::string &SEDMLFilename, const std::string &SBMLFilename, unsigned int sedmlLevel=1, unsigned int sedmlVersion=1, bool overwrite=false)
 
const std::string exportModelAndTasksToString (CCopasiDataModel &dataModel, std::string &sbmldocument, unsigned int sedmlLevel, unsigned int sedmlVersion)
 
SedDocument * getSEDMLDocument ()
 
 ~CSEDMLExporter ()
 

Protected Attributes

SedDocument * mpSEDMLDocument
 
SedUniformTimeCourse * mpTimecourse
 
SedTask * mpTimecourseTask
 
unsigned int mSEDMLLevel
 
unsigned int mSEDMLVersion
 

Detailed Description

Definition at line 37 of file CSEDMLExporter.h.

Constructor & Destructor Documentation

CSEDMLExporter::CSEDMLExporter ( )

Definition at line 649 of file CSEDMLExporter.cpp.

650  : mpSEDMLDocument(NULL)
651  , mpTimecourse(NULL)
652  , mpTimecourseTask(NULL)
653 {
654  // TODO Auto-generated constructor stub
655 }
SedDocument * mpSEDMLDocument
SedTask * mpTimecourseTask
SedUniformTimeCourse * mpTimecourse
CSEDMLExporter::~CSEDMLExporter ( )

Definition at line 657 of file CSEDMLExporter.cpp.

658 {
659  // TODO Auto-generated destructor stub
660 }

Member Function Documentation

void CSEDMLExporter::createDataGenerators ( CCopasiDataModel dataModel,
std::string &  taskId,
CCopasiTask task = NULL 
)

Creates the data generators for SEDML.

Definition at line 431 of file CSEDMLExporter.cpp.

References createDataGenerator(), CCopasiMessage::ERROR, CReportDefinition::getBodyAddr(), CPlotItem::getChannels(), CCopasiObject::getCN(), CCopasiDataModel::getDataObject(), CReportDefinition::getHeaderAddr(), CPlotSpecification::getItems(), CCopasiDataModel::getModel(), SEDMLUtils::getNextId(), CCopasiContainer::getObject(), CCopasiObject::getObjectDisplayName(), CCopasiObject::getObjectName(), CCopasiDataModel::getPlotDefinitionList(), CCopasiTask::getReport(), CReport::getReportDefinition(), CReportDefinition::getSeparator(), CReportDefinition::getTableAddr(), CReport::getTarget(), SEDMLUtils::getXPathAndName(), CPlotSpecification::isLogX(), CPlotSpecification::isLogY(), CReportDefinition::isTable(), mpSEDMLDocument, SEDMLUtils::removeCharactersFromString(), SEDML_TIME_URN, CCopasiVector< T >::size(), and CCopasiMessage::WARNING.

Referenced by createTasks().

434 {
435  const CModel* pModel = dataModel.getModel();
436  std::vector<std::string> stringsContainer; //split string container
437 
438  if (pModel == NULL)
439  CCopasiMessage(CCopasiMessage::ERROR, "SED-ML: No model for this SED-ML document. An SBML model must exist for every SED-ML document.");
440 
441  SedPlot2D* pPSedPlot;
442  SedCurve* pCurve; // = pPSedPlot->createCurve();
443 
444  //create generator for special varibale time
445  const CCopasiObject* pTime = static_cast<const CCopasiObject *>(dataModel.getModel()->getObject(CCopasiObjectName("Reference=Time")));
446  SedDataGenerator *pTimeDGenp = this->mpSEDMLDocument->createDataGenerator();
447  pTimeDGenp->setId("time");
448  pTimeDGenp->setName(pTime->getObjectName());
449  SedVariable *pTimeVar = pTimeDGenp->createVariable();
450  pTimeVar->setId("var_time");
451  pTimeVar->setTaskReference(taskId);
452  pTimeVar->setSymbol(SEDML_TIME_URN);
453  pTimeDGenp->setMath(SBML_parseFormula(pTimeVar->getId().c_str()));
454 
455  size_t i, imax = dataModel.getPlotDefinitionList()->size();
456  SedDataGenerator *pPDGen;
457 
458  if (imax == 0 && (task == NULL || task->getReport().getTarget().empty()))
459  CCopasiMessage(CCopasiMessage::ERROR, "SED-ML: No plot/report definition for this SED-ML document.");
460 
461  // export report
462  if (task != NULL && !task->getReport().getTarget().empty())
463  {
465 
466  if (def != NULL)
467  {
468  SedReport* pReport = mpSEDMLDocument->createReport();
469  std::string name = def->getObjectName();
471  //
472  pReport->setId(SEDMLUtils::getNextId("report", mpSEDMLDocument->getNumOutputs()));
473  pReport->setName(name);
474 
475  std::vector<CRegisteredObjectName> header = *def->getHeaderAddr();
476  std::vector<CRegisteredObjectName> body =
477  def->isTable() ? *def->getTableAddr() :
478  *def->getBodyAddr();
479 
480  int dsCount = 0;
481 
482  for (size_t i = 0; i < body.size(); ++i)
483  {
484  CRegisteredObjectName& current = body[i];
485 
486  if (current == def->getSeparator().getCN()) continue;
487 
488  CCopasiObject *object = dataModel.getDataObject(current);
489 
490  if (object == NULL) continue;
491 
492  const std::string& typeX = object->getObjectName();
493  std::string xAxis = object->getObjectDisplayName();
494 
495  std::string targetXPathStringX = SEDMLUtils::getXPathAndName(xAxis, typeX,
496  pModel, dataModel);
497 
498  if (object->getCN() == pTime->getCN())
499  pPDGen = pTimeDGenp;
500  else
501  pPDGen = createDataGenerator(
502  this->mpSEDMLDocument,
503  xAxis,
504  targetXPathStringX,
505  taskId,
506  i,
507  0
508  );
509 
510  SedDataSet* pDS = pReport->createDataSet();
511  pDS->setId(SEDMLUtils::getNextId("ds", ++dsCount));
512 
513  if (def->isTable())
514  {
515  CCopasiObject *headerObj = NULL;
516 
517  if (header.size() > i)
518  headerObj = dataModel.getDataObject(header[i]);
519  else
520  headerObj = dataModel.getDataObject(body[i]);
521 
522  if (headerObj != NULL)
523  pDS->setLabel(headerObj->getObjectDisplayName());
524  else
525  pDS->setLabel(xAxis);
526  }
527  else
528  pDS->setLabel(xAxis);
529 
530  pDS->setDataReference(pPDGen->getId());
531  }
532  }
533  }
534 
535  // export plots
536  for (i = 0; i < imax; i++)
537  {
538  pPSedPlot = this->mpSEDMLDocument->createPlot2D();
539  const CPlotSpecification* pPlot = (*dataModel.getPlotDefinitionList())[i];
540  std::string plotName = pPlot->getObjectName();
541 
543 
544  pPSedPlot->setId(SEDMLUtils::getNextId("plot", mpSEDMLDocument->getNumOutputs()));
545  pPSedPlot->setName(plotName);
546 
547  size_t j, jmax = pPlot->getItems().size();
548 
549  for (j = 0; j < jmax; j++)
550  {
551  const CPlotItem* pPlotItem = pPlot->getItems()[j];
552 
553  CCopasiObject *objectX, *objectY;
554 
555  if (pPlotItem->getChannels().size() >= 1)
556  {
557  objectX = dataModel.getDataObject(pPlotItem->getChannels()[0]);
558  }
559  else
560  {
561  CCopasiMessage(CCopasiMessage::WARNING, "SED-ML: Can't export plotItem '%s', as it has no data channel.", pPlotItem->getObjectName().c_str());
562  continue;
563  }
564 
565  if (objectX == NULL)
566  {
567  CCopasiMessage(CCopasiMessage::WARNING, "SED-ML: Can't export plotItem '%s' variable '%s', as it cannot be resolved.", pPlotItem->getObjectName().c_str(), pPlotItem->getChannels()[0].c_str());
568  continue;
569  }
570 
571  bool xIsTime = objectX->getCN() == pTime->getCN();
572 
573  if (pPlotItem->getChannels().size() >= 2)
574  {
575  objectY = dataModel.getDataObject(pPlotItem->getChannels()[1]);
576  }
577  else
578  {
579  CCopasiMessage(CCopasiMessage::WARNING, "SED-ML: Can't export plotItem '%s', as it has only 1 data channel.", pPlotItem->getObjectName().c_str());
580  continue;
581  }
582 
583  if (objectY == NULL)
584  {
585  CCopasiMessage(CCopasiMessage::WARNING, "SED-ML: Can't export plotItem '%s' variable '%s', as it cannot be resolved.", pPlotItem->getObjectName().c_str(), pPlotItem->getChannels()[1].c_str());
586  continue;
587  }
588 
589  const std::string& type = objectY->getObjectName();
590  std::string yAxis = objectY->getObjectDisplayName();
591  std::string sbmlId = yAxis;
592  std::string targetXPathString = SEDMLUtils::getXPathAndName(sbmlId, type,
593  pModel, dataModel);
594 
595  if (targetXPathString.empty())
596  {
597  CCopasiMessage(CCopasiMessage::WARNING, "SED-ML: Can't export plotItem '%s' variable '%s', as no xpath expression for it could be generated.", pPlotItem->getObjectName().c_str(), pPlotItem->getChannels()[1].c_str());
598  continue;
599  }
600 
601  pPDGen = createDataGenerator(
602  this->mpSEDMLDocument,
603  sbmlId,
604  targetXPathString,
605  taskId,
606  i,
607  j
608  );
609 
610  pPDGen->setName(yAxis);
611 
612  pCurve = pPSedPlot->createCurve();
613  std::ostringstream idCurveStrStream;
614  idCurveStrStream << "p";
615  idCurveStrStream << i + 1;
616  idCurveStrStream << "_curve_";
617  idCurveStrStream << j + 1;
618  pCurve->setId(idCurveStrStream.str());
619  pCurve->setLogX(pPlot->isLogX());
620  pCurve->setLogY(pPlot->isLogY());
621  pCurve->setName(yAxis);
622  pCurve->setYDataReference(pPDGen->getId());
623 
624  if (xIsTime)
625  {
626  pCurve->setXDataReference(pTimeDGenp->getId());
627  }
628  else
629  {
630  const std::string& typeX = objectX->getObjectName();
631  std::string xAxis = objectX->getObjectDisplayName();
632  std::string targetXPathStringX = SEDMLUtils::getXPathAndName(xAxis, typeX,
633  pModel, dataModel);
634 
635  pPDGen = createDataGenerator(
636  this->mpSEDMLDocument,
637  xAxis,
638  targetXPathStringX,
639  taskId,
640  i,
641  j
642  );
643  pCurve->setXDataReference(pPDGen->getId());
644  }
645  }
646  }
647 }
CCopasiObject * getDataObject(const CCopasiObjectName &CN) const
virtual std::string getObjectDisplayName(bool regular=true, bool richtext=false) const
virtual CCopasiObjectName getCN() const
const std::string & getObjectName() const
virtual size_t size() const
SedDataGenerator * createDataGenerator(SedDocument *mpSEDMLDocument, const std::string &sbmlId, const std::string &targetXPathString, const std::string &taskId, size_t i, size_t j)
SedDocument * mpSEDMLDocument
const CCopasiReportSeparator & getSeparator() const
std::vector< CPlotDataChannelSpec > & getChannels()
Definition: CPlotItem.cpp:214
const std::string & getTarget() const
Definition: CReport.cpp:89
static std::string getXPathAndName(std::string &sbmlId, const std::string &type, const CModel *pModel, const CCopasiDataModel &dataModel)
Definition: SEDMLUtils.cpp:65
std::vector< CRegisteredObjectName > * getTableAddr()
const COutputDefinitionVector * getPlotDefinitionList() const
CReport & getReport()
const CCopasiVector< CPlotItem > & getItems() const
std::vector< CRegisteredObjectName > * getBodyAddr()
static std::string getNextId(const std::string &base, int count)
Definition: SEDMLUtils.cpp:294
Definition: CModel.h:50
virtual const CObjectInterface * getObject(const CCopasiObjectName &cn) const
CReportDefinition * getReportDefinition()
Definition: CReport.cpp:83
static std::string & removeCharactersFromString(std::string &str, const std::string &characters)
Definition: SEDMLUtils.cpp:273
std::vector< CRegisteredObjectName > * getHeaderAddr()
#define SEDML_TIME_URN
Definition: SEDMLUtils.h:21
void CSEDMLExporter::createModels ( CCopasiDataModel dataModel,
std::string &  modelRef 
)

Creates the models for SEDML.

Definition at line 362 of file CSEDMLExporter.cpp.

References mpSEDMLDocument.

Referenced by createSEDMLDocument().

363 {
364  SedModel *model = this->mpSEDMLDocument->createModel();
365  model->setId(modelRef.substr(0, modelRef.length() - 4));
366  model->setSource(modelRef);
367  model->setLanguage("urn:sedml:language:sbml");
368 }
SedDocument * mpSEDMLDocument
std::string CSEDMLExporter::createScanTask ( CCopasiDataModel dataModel,
const std::string &  modelId 
)

Creates a scan task for the given model id if the dataModel contains a number of scan items and returns its id.

Returns
the id of the task, if craeted, otherwise an empty string.

Creates the simulations for SEDML.

Definition at line 181 of file CSEDMLExporter.cpp.

References createSteadyStateTask(), CScanProblem::getContinueFromCurrentState(), SEDMLUtils::getNextId(), CScanProblem::getNumberOfScanItems(), CCopasiContainer::getObject(), CCopasiParameterGroup::getParameter(), CCopasiTask::getProblem(), CScanProblem::getScanItem(), CScanProblem::getSubtask(), CCopasiDataModel::getTaskList(), CCopasiParameter::getValue(), SEDMLUtils::getXPathForObject(), max, min, mpSEDMLDocument, mpTimecourseTask, CCopasiParameter::Value::pBOOL, CCopasiParameter::Value::pCN, CCopasiParameter::Value::pDOUBLE, pTask, CCopasiParameter::Value::pUINT, CScanProblem::SCAN_LINEAR, CScanProblem::SCAN_RANDOM, CScanProblem::SCAN_REPEAT, SEDML_TIME_URN, CCopasiTask::steadyState, CCopasiTask::timeCourse, and CCopasiMessage::WARNING.

Referenced by createTasks().

182 {
183  // need L1V2 to export repeated tasks
184  if (mpSEDMLDocument->getVersion() != 2) return "";
185 
186  CScanTask* pTask = dynamic_cast<CScanTask*>((*dataModel.getTaskList())["Scan"]);
187 
188  if (pTask == NULL) return "";
189 
190  CScanProblem* pProblem = dynamic_cast<CScanProblem*>(pTask->getProblem());
191  size_t numItems = pProblem->getNumberOfScanItems();
192 
193  if (numItems == 0)
194  return "";
195 
196  if (pProblem->getSubtask() != CCopasiTask::steadyState &&
197  pProblem->getSubtask() != CCopasiTask::timeCourse)
198  {
199  CCopasiMessage(CCopasiMessage::WARNING, "SED-ML: This version of COPASI only allows the export of time course or steady state scans.");
200  return "";
201  }
202 
203  std::string subTaskId;
204 
205  if (pProblem->getSubtask() == CCopasiTask::steadyState)
206  {
207  subTaskId = createSteadyStateTask(dataModel, modelId);
208  }
209  else
210  {
211  subTaskId = mpTimecourseTask->getId();
212  }
213 
214  SedRepeatedTask* task = mpSEDMLDocument->createRepeatedTask();
215  std::string taskId = SEDMLUtils::getNextId("task", mpSEDMLDocument->getNumTasks());
216  task->setId(taskId);
217  task->setResetModel(!pProblem->getContinueFromCurrentState());
218 
219  // craete ranges / changes
220  for (size_t i = 0; i < numItems; ++i)
221  {
222  CCopasiParameterGroup* current = pProblem->getScanItem(i);
223  CScanProblem::Type type = (CScanProblem::Type)(*current->getParameter("Type")->getValue().pUINT);
224 
225  // ignore random items
226  if (type == CScanProblem::SCAN_RANDOM)
227  {
228  CCopasiMessage(CCopasiMessage::WARNING, "SED-ML: This version of COPASI cannot export random scan items, they will be ignored.");
229  continue;
230  }
231 
232  int numSteps = (*current->getParameter("Number of steps")->getValue().pUINT);
233 
234  // handle repeats
235  if (type == CScanProblem::SCAN_REPEAT)
236  {
237  SedUniformRange *range = task->createUniformRange();
238  range->setId(SEDMLUtils::getNextId("range", task->getNumRanges()));
239  range->setStart(0);
240  range->setEnd(numSteps);
241  range->setNumberOfPoints(numSteps);
242  range->setType("linear");
243 
244  if (task->isSetRangeId())
245  task->setRangeId(range->getId());
246 
247  continue;
248  }
249 
250  // handle scans
251  if (type == CScanProblem::SCAN_LINEAR)
252  {
253  double min = (*current->getParameter("Minimum")->getValue().pDOUBLE);
254  double max = (*current->getParameter("Maximum")->getValue().pDOUBLE);
255  bool log = (*current->getParameter("log")->getValue().pBOOL);
256 
257  SedUniformRange *range = task->createUniformRange();
258  range->setId(SEDMLUtils::getNextId("range", task->getNumRanges()));
259  range->setStart(min);
260  range->setEnd(max);
261  range->setNumberOfPoints(numSteps);
262  range->setType(log ? "log" : "linear");
263 
264  const CRegisteredObjectName& cn = (*current->getParameter("Object")->getValue().pCN);
265  std::string xpath = SEDMLUtils::getXPathForObject(*static_cast<const CCopasiObject*>(dataModel.getObject(cn)));
266 
267  if (xpath.empty())
268  {
269  CCopasiMessage(CCopasiMessage::WARNING, "SED-ML: This version of COPASI cannot export the selected scan object, it will be ignored.");
270  continue;
271  }
272 
273  SedSetValue *change = task->createTaskChange();
274  change->setModelReference(modelId);
275 
276  if (xpath == SEDML_TIME_URN)
277  {
278  change->setSymbol(xpath);
279  }
280  else
281  {
282  change->setTarget(xpath);
283  }
284 
285  change->setRange(range->getId());
286  change->setMath(SBML_parseFormula(range->getId().c_str()));
287 
288  continue;
289  }
290  }
291 
292  if (!task->isSetRangeId() && task->getNumRanges() > 0)
293  task->setRangeId(task->getRange(0)->getId());
294 
295  // create subtask
296  SedSubTask* subTask = task->createSubTask();
297  subTask->setOrder(1);
298  subTask->setTask(subTaskId);
299 
300  return taskId;
301 }
CCopasiProblem * getProblem()
SedDocument * mpSEDMLDocument
size_t getNumberOfScanItems() const
CRegisteredObjectName * pCN
CTSSATask * pTask
SedTask * mpTimecourseTask
const Value & getValue() const
CCopasiVectorN< CCopasiTask > * getTaskList()
unsigned C_INT32 * pUINT
CCopasiParameter * getParameter(const std::string &name)
bool getContinueFromCurrentState() const
static std::string getXPathForObject(const CCopasiObject &object)
Definition: SEDMLUtils.cpp:284
static std::string getNextId(const std::string &base, int count)
Definition: SEDMLUtils.cpp:294
CCopasiTask::Type getSubtask() const
const CCopasiParameterGroup * getScanItem(size_t index) const
virtual const CObjectInterface * getObject(const CCopasiObjectName &cn) const
std::string createSteadyStateTask(CCopasiDataModel &dataModel, const std::string &modelId)
#define SEDML_TIME_URN
Definition: SEDMLUtils.h:21
#define min(a, b)
Definition: f2c.h:175
#define max(a, b)
Definition: f2c.h:176
void CSEDMLExporter::createSEDMLDocument ( CCopasiDataModel dataModel,
std::string  modelRef 
)

Creates an SEDMLDocument and SBMLDocument from the given CCopasiDataModelObject. It checks if an SEDMLDocument and SBMLDocument already exists from an import and if this is the case, the old document is copied. If none exists a new one is created. Copying the old one makes sure that if something goes wrong during export, the original model is still consistent.

Definition at line 155 of file CSEDMLExporter.cpp.

References createModels(), createTasks(), fatalError, CCopasiDataModel::getModel(), CCopasiDataModel::getPlotDefinitionList(), mpSEDMLDocument, mSEDMLLevel, and mSEDMLVersion.

Referenced by exportModelAndTasksToString().

156 {
157  const SedDocument* pOldSEDMLDocument = NULL; //dataModel.getCurrentSEDMLDocument();
158  const CModel* pModel = dataModel.getModel();
159  COutputDefinitionVector *plotDef = dataModel.getPlotDefinitionList();
160  assert(plotDef != NULL); //need to emit a message
161  assert(pModel != NULL);
162 
163  if (pOldSEDMLDocument == NULL)
164  {
165  this->mpSEDMLDocument = new SedDocument(mSEDMLLevel, mSEDMLVersion);
166  }
167  else
168  {
169  this->mpSEDMLDocument = dynamic_cast<SedDocument*>(pOldSEDMLDocument->clone());
170  }
171 
172  if (this->mpSEDMLDocument == NULL) fatalError();
173 
174  createModels(dataModel, modelRef);
175  createTasks(dataModel, modelRef);
176 }
unsigned int mSEDMLVersion
SedDocument * mpSEDMLDocument
#define fatalError()
unsigned int mSEDMLLevel
void createTasks(CCopasiDataModel &dataModel, std::string &modelRef)
void createModels(CCopasiDataModel &dataModel, std::string &modelRef)
const COutputDefinitionVector * getPlotDefinitionList() const
Definition: CModel.h:50
std::string CSEDMLExporter::createSteadyStateTask ( CCopasiDataModel dataModel,
const std::string &  modelId 
)

Creates the steady state task for the given model id and returns its id.

Creates the simulations for SEDML.

Definition at line 338 of file CSEDMLExporter.cpp.

References SEDMLUtils::getNextId(), CCopasiTask::getProblem(), CCopasiDataModel::getTaskList(), mpSEDMLDocument, and pTask.

Referenced by createScanTask().

339 {
340  SedSteadyState *steady = this->mpSEDMLDocument->createSteadyState();
341  steady->setId(SEDMLUtils::getNextId("steady", mpSEDMLDocument->getNumSimulations()));
342  //presently SEDML only supports time course
343  CCopasiTask* pTask = (*dataModel.getTaskList())["Steady-State"];
344  CTrajectoryProblem* tProblem = static_cast<CTrajectoryProblem*>(pTask->getProblem());
345 
346  // set the correct KISAO Term
347  SedAlgorithm* alg = steady->createAlgorithm();
348  alg->setKisaoID("KISAO:0000282");
349 
350  SedTask *task = this->mpSEDMLDocument->createTask();
351  std::string taskId = SEDMLUtils::getNextId("task", mpSEDMLDocument->getNumTasks());
352  task->setId(taskId);
353  task->setSimulationReference(steady->getId());
354  task->setModelReference(modelId);
355 
356  return taskId;
357 }
CCopasiProblem * getProblem()
SedDocument * mpSEDMLDocument
CTSSATask * pTask
CCopasiVectorN< CCopasiTask > * getTaskList()
static std::string getNextId(const std::string &base, int count)
Definition: SEDMLUtils.cpp:294
void CSEDMLExporter::createTasks ( CCopasiDataModel dataModel,
std::string &  modelRef 
)

Creates the Tasks for SEDML.

Creates the Tasks for SEDML. This will always create a task running a time course simulation. If the parameter scan has been specified, it will be exported as well.

Definition at line 374 of file CSEDMLExporter.cpp.

References createDataGenerators(), createScanTask(), createTimeCourseTask(), and CCopasiDataModel::getTaskList().

Referenced by createSEDMLDocument().

375 {
376  std::string modelId = modelRef.substr(0, modelRef.length() - 4);
377  // create time course task
378  std::string taskId = createTimeCourseTask(dataModel, modelId);
379  createDataGenerators(dataModel, taskId, (*dataModel.getTaskList())["Time-Course"]);
380 
381  taskId = createScanTask(dataModel, modelId);
382 
383  if (!taskId.empty())
384  createDataGenerators(dataModel, taskId, (*dataModel.getTaskList())["Scan"]);
385 }
std::string createScanTask(CCopasiDataModel &dataModel, const std::string &modelId)
void createDataGenerators(CCopasiDataModel &dataModel, std::string &taskId, CCopasiTask *task=NULL)
std::string createTimeCourseTask(CCopasiDataModel &dataModel, const std::string &modelId)
CCopasiVectorN< CCopasiTask > * getTaskList()
std::string CSEDMLExporter::createTimeCourseTask ( CCopasiDataModel dataModel,
const std::string &  modelId 
)

Creates the timecourse task for the given model id and returns its id.

Creates the simulations for SEDML.

Definition at line 306 of file CSEDMLExporter.cpp.

References CCopasiTask::getMethod(), SEDMLUtils::getNextId(), CCopasiObject::getObjectName(), CTrajectoryProblem::getOutputStartTime(), CCopasiTask::getProblem(), CTrajectoryProblem::getStepNumber(), CTrajectoryProblem::getStepSize(), CCopasiDataModel::getTaskList(), mpSEDMLDocument, mpTimecourse, mpTimecourseTask, and pTask.

Referenced by createTasks().

307 {
308  mpTimecourse = this->mpSEDMLDocument->createUniformTimeCourse();
309  mpTimecourse->setId(SEDMLUtils::getNextId("sim", mpSEDMLDocument->getNumSimulations()));
310  //presently SEDML only supports time course
311  CCopasiTask* pTask = (*dataModel.getTaskList())["Time-Course"];
312  CTrajectoryProblem* tProblem = static_cast<CTrajectoryProblem*>(pTask->getProblem());
313  mpTimecourse->setInitialTime(0.0);
314  mpTimecourse->setOutputStartTime(tProblem->getOutputStartTime());
315  mpTimecourse->setOutputEndTime(tProblem->getStepNumber()*tProblem->getStepSize());
316  mpTimecourse->setNumberOfPoints(tProblem->getStepNumber());
317 
318  // set the correct KISAO Term
319  SedAlgorithm* alg = mpTimecourse->createAlgorithm();
320 
321  if (pTask->getMethod()->getObjectName().find("Stochastic") != std::string::npos)
322  alg->setKisaoID("KISAO:0000241");
323  else
324  alg->setKisaoID("KISAO:0000019");
325 
326  mpTimecourseTask = this->mpSEDMLDocument->createTask();
327  std::string taskId = SEDMLUtils::getNextId("task", mpSEDMLDocument->getNumTasks());
328  mpTimecourseTask->setId(taskId);
329  mpTimecourseTask->setSimulationReference(mpTimecourse->getId());
330  mpTimecourseTask->setModelReference(modelId);
331 
332  return taskId;
333 }
CCopasiProblem * getProblem()
const std::string & getObjectName() const
const unsigned C_INT32 & getStepNumber() const
SedDocument * mpSEDMLDocument
CTSSATask * pTask
const C_FLOAT64 & getStepSize() const
SedTask * mpTimecourseTask
CCopasiVectorN< CCopasiTask > * getTaskList()
SedUniformTimeCourse * mpTimecourse
CCopasiMethod * getMethod()
const C_FLOAT64 & getOutputStartTime() const
static std::string getNextId(const std::string &base, int count)
Definition: SEDMLUtils.cpp:294
bool CSEDMLExporter::exportModelAndTasks ( CCopasiDataModel dataModel,
const std::string &  filename,
const std::string &  sbmlDocument,
unsigned int  sedmlLevel = 1,
unsigned int  sedmlVersion = 1,
bool  overwrite = false 
)

Export the model and Task to SEDML. The SEDML document is written to the file given by SEDMLFilename and reference SBML model is written to SBMLFilename . If the export fails, false is returned.

Definition at line 92 of file CSEDMLExporter.cpp.

References createUniqueModelFileName(), CDirEntry::dirName(), CCopasiMessage::ERROR, exportModelAndTasksToString(), CLocaleString::fromUtf8(), MCDirEntry, and CDirEntry::Separator.

98 {
99  bool success = true;
100  /* create a string that represents the SBMLDocument */
101 
102  // std::string sedmlModelSource = "model1.xml"; //always name of the SBML model to reference in SEDML document
103  // create a unique name for all exported models, rather than to overwrite existing ones!
104  std::string sedmlModelSource = createUniqueModelFileName(CDirEntry::dirName(filename), "model", ".xml");
105 
106  std::string sbmlFileName;
107  sbmlFileName = CDirEntry::dirName(filename) + CDirEntry::Separator + sedmlModelSource;
108 
109  std::ifstream sbmlFile(CLocaleString::fromUtf8(sbmlFileName).c_str(), std::ios::in);
110 
111  if (sbmlFile && !overwrite)
112  {
113  // create a CCopasiMessage with the appropriate error
114  CCopasiMessage(CCopasiMessage::ERROR, MCDirEntry + 1, sbmlFileName.c_str());
115  return false;
116  }
117 
118  /* write the sbml model document to a file */
119  std::ofstream sbmlOutFile(CLocaleString::fromUtf8(sbmlFileName).c_str(), std::ios::out | std::ios::trunc);
120  sbmlOutFile << sbmlDocument;
121  sbmlOutFile.close();
122 
123  std::string str = this->exportModelAndTasksToString(dataModel, sedmlModelSource, sedmlLevel, sedmlVersion);
124 
125  //std::cout<<str<<std::endl; //only for debuging
126 
127  if (!str.empty())
128  {
129  /* check if the file already exists.
130  If yes, write if overwrite is true,
131  else create an appropriate CCopasiMessage. */
132  std::ifstream testInfile(CLocaleString::fromUtf8(filename).c_str(), std::ios::in);
133 
134  if (testInfile && !overwrite)
135  {
136  // create a CCopasiMessage with the appropriate error
137  CCopasiMessage(CCopasiMessage::ERROR, MCDirEntry + 1, filename.c_str());
138  return false;
139  }
140 
141  /* write the document to a file */
142  std::ofstream outfile(CLocaleString::fromUtf8(filename).c_str(), std::ios::out | std::ios::trunc);
143  outfile << str;
144  outfile.close();
145  }
146  else
147  {
148  /* if no SBMLDocument could be created return false */
149  success = false;
150  }
151 
152  return success;
153 }
const std::string exportModelAndTasksToString(CCopasiDataModel &dataModel, std::string &sbmldocument, unsigned int sedmlLevel, unsigned int sedmlVersion)
static std::string dirName(const std::string &path)
Definition: CDirEntry.cpp:135
std::string createUniqueModelFileName(const std::string &dir, const std::string &baseName, const std::string &extension=".xml")
static const std::string Separator
Definition: CDirEntry.h:34
#define MCDirEntry
static CLocaleString fromUtf8(const std::string &utf8)
const std::string CSEDMLExporter::exportModelAndTasksToString ( CCopasiDataModel dataModel,
std::string &  sbmlModelSource,
unsigned int  sedmlLevel,
unsigned int  sedmlVersion 
)

Export the model and Task to SEDML. The SEDML document is returned as a string and SBML model is copied to sbmldocument parameter. In case of an error, an empty string is returned.

Definition at line 46 of file CSEDMLExporter.cpp.

References createSEDMLDocument(), mpSEDMLDocument, mSEDMLLevel, mSEDMLVersion, pdelete, and CVersion::VERSION.

Referenced by exportModelAndTasks().

50 {
51  this->mSEDMLLevel = sedmlLevel;
52  this->mSEDMLVersion = sedmlVersion;
53  this->createSEDMLDocument(dataModel, sbmlModelSource);
54 
55  CSBMLExporter exporter;
56  SedWriter* writer = new SedWriter();
57 
58  writer->setProgramName("COPASI");
59  writer->setProgramVersion(CVersion::VERSION.getVersion().c_str());
60 
61  char* d = writer->writeToString(this->mpSEDMLDocument);
62  std::string returnValue = d;
63 
64  if (d) free(d);
65 
66  pdelete(writer);
67 
68  return returnValue;
69 }
unsigned int mSEDMLVersion
#define pdelete(p)
Definition: copasi.h:215
SedDocument * mpSEDMLDocument
unsigned int mSEDMLLevel
void createSEDMLDocument(CCopasiDataModel &dataModel, std::string modelRef)
static const CVersion VERSION
Definition: CVersion.h:164
SedDocument* CSEDMLExporter::getSEDMLDocument ( )
inline

Definition at line 51 of file CSEDMLExporter.h.

References mpSEDMLDocument.

51 {return this->mpSEDMLDocument;};
SedDocument * mpSEDMLDocument

Member Data Documentation

SedDocument* CSEDMLExporter::mpSEDMLDocument
protected
SedUniformTimeCourse* CSEDMLExporter::mpTimecourse
protected

Definition at line 44 of file CSEDMLExporter.h.

Referenced by createTimeCourseTask().

SedTask* CSEDMLExporter::mpTimecourseTask
protected

Definition at line 45 of file CSEDMLExporter.h.

Referenced by createScanTask(), and createTimeCourseTask().

unsigned int CSEDMLExporter::mSEDMLLevel
protected

Definition at line 42 of file CSEDMLExporter.h.

Referenced by createSEDMLDocument(), and exportModelAndTasksToString().

unsigned int CSEDMLExporter::mSEDMLVersion
protected

Definition at line 43 of file CSEDMLExporter.h.

Referenced by createSEDMLDocument(), and exportModelAndTasksToString().


The documentation for this class was generated from the following files: