## Compucell3D files for
## "3D multi-cell simulation of tumor growth and angiogenesis."
## Shirinifard A, Gens JS, Zaitlen BL, Poplawski NJ, Swat M, Glazier JA.
## PLoS One. 2009 Oct 16;4(10):e7190.
## 
## http://www.ncbi.nlm.nih.gov/pubmed/19834621
## http://www.ncbi.nlm.nih.gov/pmc/articles/PMC2760204/
## 
## Tumor growth and angiogenesis computational model
## Authors:
##   Shirinifard A, Gens JS, Zaitlen BL, Poplawski NJ, Swat M, Glazier JA.
##   Indiana University, Bloomington
## Computational Platform:
##   CompuCell3D CC3D
##   
## Foundational Model of Anatomy FMA:
##   GO:0001525 angiogenesis
##   GO:0045766 activation of angiogenesis
## 
## Gene Ontology GO:
##   GO:0010573 VEGF production 
##   GO:0051301 cell division
##   GO:0008219 cell death necrosis
##   GO:0016049 cell growth
##   GO:0071456 "cellular response to hypoxia" "cellular response to hypoxic stress"
##   GO:0001666 "response to hypoxia"
##   GO:0007067 mitosis
##   
## Chemical Entities of Biological Interest ChEBI: 
##   CHEBI:15379 oxygen dioxygen
## 
## Molecular Interactions MI:
##   MI:2111 "diffusion coeff" "diffusion coefficient"
## 
## Protein Ontology PRO:
##   PRO:000017284 	VEGF-A 
## 
## Mammalian Phenotype MA:
##   MP:0005039 hypoxia
##   
## EMBRACE Ontology for Data and Methods EDAM:
##   EDAM:0000950 Biological model
## 
## Medical Subject Heading MeSH:
##   MeSH:C04:Neoplasms neoplasm
##   MeSH:C23.550.589.500 Neovascularization, Pathologic
##   MeSH:G09.330.190.751 Neovascularization, Physiologic
##   MeSH:E05.599 Models, Theoretical
##   MeSH:E05.599.395 Models, Biological
##   MeSH:L01.224.160 Computer Simulation
## 
## Computational Platform:
##   CompuCell3D CC3D
##   cellular potts model
## 
## HERBS AND JIMS BIG DATA AND COMPUTATIONAL MODEL REPOSITORY
## Indiana University Bloomington
## University of Washington
## 
## Three python scripts, delineated by ## ++++++..., for CC3D are included below;
##  angio_growth_08052009_01_45_36.py
##  angio_growth_plugins_08052009_01_45_36.py
##  angio_growth_steppables_08052009_01_45_36.py
##  
## ++++++++++++++++++++++++++++++
## angio_growth_08052009_01_45_36.py
## ++++++++++++++++++++++++++++++
####
#### The simulation code is compatible with CompuCell3D ver 3.3.1
####

import sys
from os import environ
import string
sys.path.append(environ["PYTHON_MODULE_PATH"])

#

import CompuCellSetup

sim,simthread = CompuCellSetup.getCoreSimulationObjects()

#Create extra player fields here or add attributes

pyAttributeAdder,listAdder=CompuCellSetup.attachListToCells(sim)

CompuCellSetup.initializeSimulationObjects(sim,simthread)

#
##################
##########	PLUGINS
##################
#

import CompuCell

from angio_growth_plugins_08052009_01_45_36 import *

changeWatcherRegistry=CompuCellSetup.getChangeWatcherRegistry(sim)

stepperRegistry=CompuCellSetup.getStepperRegistry(sim)

mitPy=MitosisPyPlugin(sim,changeWatcherRegistry,stepperRegistry)

#### seting doubling volumes for normal, hypoxic, ActiveNeovascular, InactiveNeovascular
doublingVolumeDict = {1:54,2:54,4:80,6:80}
mitPy.setCellDoublingVolume(doublingVolumeDict)

#
##################
##########	STEPPABLES
##################
#

from PySteppables import SteppableRegistry
steppableRegistry=SteppableRegistry()

from angio_growth_steppables_08052009_01_45_36 import *
                                         #sim,frequency,areaThresh,nutrientThresh,necroticThresh
volumeParamSteppable=VolumeParamSteppable(sim,1,1,5,1)
steppableRegistry.registerSteppable(volumeParamSteppable)

#
##################
##########	COMPUCELL3D LOOPS
##################
#

CompuCellSetup.mainLoop(sim,simthread,steppableRegistry)

## ++++++++++++++++++++++++++++++
## angio_growth_plugins_08052009_01_45_36.py
## ++++++++++++++++++++++++++++++
####
#### The simulation code is compatible with CompuCell3D ver 3.3.1
####

from CompuCell import MitosisSimplePlugin      
from PyPlugins import *
from PySteppables import CellList
from CompuCell import NeighborFinderParams
import time,sys

class MitosisPyPluginBase(StepperPy,Field3DChangeWatcherPy):
   def __init__(self,_simulator,_changeWatcherRegistry,_stepperRegistry):

      Field3DChangeWatcherPy.__init__(self,_changeWatcherRegistry)
      self.simulator=_simulator
      self.mitosisPlugin=MitosisSimplePlugin()
      self.mitosisPlugin.setPotts(self.simulator.getPotts())
      self.mitosisPlugin.turnOn()
      self.mitosisPlugin.init(self.changeWatcher.sim)
      self.counter=0
      self.mitosisFlag=0
      self.doublingVolumeDict=0
      _changeWatcherRegistry.registerPyChangeWatcher(self)
      _stepperRegistry.registerPyStepper(self)        
      
   def setPotts(self,potts):
      self.mitosisPlugin.setPotts(potts)
   
   def setDoublingVolume(self,_doublingVolume):
      self.doublingVolume=_doublingVolume;
      self.mitosisPlugin.setDoublingVolume(self.doublingVolume)
   
   def setCellDoublingVolume(self,_doublingVolumeDict):
      self.doublingVolumeDict=_doublingVolumeDict;
      for i in self.doublingVolumeDict.keys():
         print self.doublingVolumeDict[i]

   def field3DChange(self):
      cell = self.changeWatcher.newCell
      if cell and self.doublingVolumeDict.has_key(cell.type) and cell.volume>self.doublingVolumeDict[cell.type]:
         print "Type: ", cell.type, " Doubling Volume: ", self.doublingVolumeDict[cell.type], " Current Volume: ", cell.volume
         self.setDoublingVolume(self.doublingVolumeDict[cell.type])
         self.mitosisPlugin.field3DChange(self.changeWatcher.changePoint,self.changeWatcher.newCell,self.changeWatcher.newCell)
         self.mitosisFlag=1
         
   def step(self):
      if self.mitosisFlag:
         print "ABOUT TO DO MITOSIS"
         self.mitosisFlag=self.mitosisPlugin.doMitosis()
         self.childCell=self.mitosisPlugin.getChildCell()
         self.parentCell=self.mitosisPlugin.getParentCell()
         self.updateAttributes()
         self.mitosisFlag=0
         
   def updateAttributes(self):
      self.childCell.targetVolume=self.parentCell.targetVolume
      self.childCell.lambdaVolume=self.parentCell.lambdaVolume
      self.childCell.type=self.parentCell.type


class MitosisPyPlugin(MitosisPyPluginBase):
   def __init__(self , _simulator , _changeWatcherRegistry , _stepperRegistry):
      MitosisPyPluginBase.__init__(self,_simulator,_changeWatcherRegistry,_stepperRegistry)

   def updateAttributes(self):
## Mitosis of normal tumor and hypoxic cells    
      if self.parentCell.type==1 or self.parentCell.type==2:
         self.childCell.type=1
	 self.childCell.targetVolume=33
         self.childCell.lambdaVolume=10
         self.childCell.targetSurface=90
         self.childCell.lambdaSurface=2
         self.parentCell.targetVolume=33
         self.parentCell.lambdaVolume=10
         self.parentCell.targetSurface=90
         self.parentCell.lambdaSurface=2
## Mitosis of ActiveNeovascular and InactiveNeovascular cells
      if self.parentCell.type==6 or self.parentCell.type==4:
         self.childCell.type=4
	 self.childCell.targetVolume=60
         self.childCell.lambdaVolume=13
         self.childCell.targetSurface=150
         self.childCell.lambdaSurface=3
         self.parentCell.targetVolume=60
         self.parentCell.lambdaVolume=13
         self.parentCell.targetSurface=150
         self.parentCell.lambdaSurface=3

## +++++++++++++++++++++++++
## angio_growth_steppables_08052009_01_45_36.py
## +++++++++++++++++++++++++
####
#### The simulation code is compatible with CompuCell3D ver 3.3.1
####
from PySteppables import *
import CompuCell
import sys
import time

class VolumeParamSteppable(SteppablePy):
   def __init__(self,_simulator,_frequency=1,_areaThresh=0,_nutrientThresh=0,_necroticThresh=0):
      SteppablePy.__init__(self,_frequency)
      self.simulator=_simulator
      self.inventory=self.simulator.getPotts().getCellInventory()
      self.cellList=CellList(self.inventory)
      self.nTrackerPlugin=CompuCell.getNeighborTrackerPlugin()
      self.areaThresh = _areaThresh
      self.nutrientThresh = _nutrientThresh
      self.necroticThresh = _necroticThresh
      self.fieldNameNeoVascular = 'VEGF2'
      self.fieldNameNormal = 'Oxygen'
      #self.output_file = open("CellDiffusionData_08052009_01_45_36.txt",'w')
      
   def start(self):
      for cell in self.cellList:
         if cell.type==4 or cell.type==5 or cell.type==6:
	    cell.targetVolume=60
            cell.lambdaVolume=13.0
	    cell.targetSurface=150
            cell.lambdaSurface=3.0
         else:
	    cell.targetVolume=33.0
            cell.lambdaVolume=10.0
	    cell.targetSurface=90.0
            cell.lambdaSurface=2
	 
   def step(self,mcs):
      fieldNeoVasc=CompuCell.getConcentrationField(self.simulator,self.fieldNameNeoVascular)
      fieldMalig=CompuCell.getConcentrationField(self.simulator,self.fieldNameNormal)
      #print mcs
      
      for cell in self.cellList:

         # Inactive neovascular differentiation
         if cell.type == 6:
            totalArea = 0
            pt=CompuCell.Point3D()
            pt.x=int(round(cell.xCM/max(float(cell.volume),0.001)))
            pt.y=int(round(cell.yCM/max(float(cell.volume),0.001)))
            pt.z=int(round(cell.zCM/max(float(cell.volume),0.001)))
            concentration=fieldNeoVasc.get(pt)
            if concentration>0.5:
            
	      cellNeighborList=CellNeighborListAuto(self.nTrackerPlugin,cell)
 	      for neighborSurfaceData in cellNeighborList:
		  #Check to ensure cell neighbor is not medium
	          if neighborSurfaceData.neighborAddress:
	 	     if neighborSurfaceData.neighborAddress.type == 5 or neighborSurfaceData.neighborAddress.type == 6 or neighborSurfaceData.neighborAddress.type == 7:
			   
			#sum up common surface area of cell with its neighbors
		        totalArea+=neighborSurfaceData.commonSurfaceArea 
			#print "concentration: ", concentration,"  commonSurfaceArea:",neighborSurfaceData.commonSurfaceArea
	      print cell.type,totalArea     
	      if totalArea < 70:
		 #Growth rate equation
		 #print cell.type,"##surface area",cell.surface,"##cell volume:",cell.volume,"##cell target volume:",cell.targetVolume,"##common surface area:",totalArea
		 cell.targetVolume+=0.06*concentration/(0.5 + concentration)
		 cell.targetSurface+=0.15*concentration/(0.5 + concentration)
	         #print 0.02*concentration/(0.5 + concentration)+0.04
         
         ## Active neovascular growth
         if cell.type == 4:
            totalArea = 0
            pt=CompuCell.Point3D()
            pt.x=int(round(cell.xCM/max(float(cell.volume),0.00000001)))
            pt.y=int(round(cell.yCM/max(float(cell.volume),0.00000001)))
            pt.z=int(round(cell.zCM/max(float(cell.volume),0.00000001)))
            concentration=fieldNeoVasc.get(pt)
            if concentration>0.5:
            
	      cellNeighborList=CellNeighborListAuto(self.nTrackerPlugin,cell)
 	      for neighborSurfaceData in cellNeighborList:
		  #Check to ensure cell neighbor is not medium
	          if neighborSurfaceData.neighborAddress:
	 	     if neighborSurfaceData.neighborAddress.type == 5 or neighborSurfaceData.neighborAddress.type == 7 or neighborSurfaceData.neighborAddress.type == 6:
			   
			#sum up common surface area of cell with its neighbors
		        totalArea+=neighborSurfaceData.commonSurfaceArea 
			#print "concentration: ", concentration,"  commonSurfaceArea:",neighborSurfaceData.commonSurfaceArea
	      #print cell.type,totalArea      
	      if totalArea < 50:
		 #Growth rate equation
		 #print cell.type,"##surface area",cell.surface,"##cell volume:",cell.volume,"##cell target volume:",cell.targetVolume,"##common surface area:",totalArea
		 cell.targetVolume+=0.06*concentration/(0.5 + concentration)
		 cell.targetSurface+=0.15*concentration/(0.5 + concentration)
	         ##print 0.02*concentration/(0.5 + concentration)+0.04

	 
	    
         #Malignat and Hypoxic Cells growth
         if cell.type == 1 or cell.type == 2:
            #print cell.volume

            pt=CompuCell.Point3D()
            pt.x=int(round(cell.xCM/max(float(cell.volume),0.001)))
            pt.y=int(round(cell.yCM/max(float(cell.volume),0.001)))
            pt.z=int(round(cell.zCM/max(float(cell.volume),0.001)))
            #self.output_file.write("%f %f %f " %(cell.xCM/cell.volume, cell.yCM/cell.volume,cell.zCM/cell.volume))

            concentration2=fieldMalig.get(pt)
            #switch to Hypoxic cell type
            if (concentration2 < self.nutrientThresh and mcs>100):
	       cell.type=2

	    #switch to Necrotic cell type
	    if (concentration2 < self.necroticThresh and mcs>100):
               cell.type=3
               
            #set growth rate equation
            if (mcs>100):
	       cell.targetVolume+=0.04*concentration2/(10+concentration2)
               cell.targetSurface+=0.12*concentration2/(10+concentration2)

         #Hypoxic Cells
         if cell.type == 2:
            #print " #Hypoxic Volume: ", cell.volume
            pt=CompuCell.Point3D()
            pt.x=int(round(cell.xCM/max(float(cell.volume),0.001)))
            pt.y=int(round(cell.yCM/max(float(cell.volume),0.001)))
            pt.z=int(round(cell.zCM/max(float(cell.volume),0.001)))
            concentration3=fieldMalig.get(pt)
            #switch to Necrotic cell type
            if (concentration3 < self.necroticThresh and mcs>100):
               cell.type=3
	    #switch to Normal cell type
	    if (concentration3 > self.nutrientThresh):
               cell.type=1
	    
         #Necrotic Cells
         if cell.type == 3:
            #set growth rate equation
            cell.targetVolume-=0.5
            cell.lambdaSurface=0
	 
      #self.output_file.write("\n")