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

# CBO CBO_Object SystemQuality ComputationalPlatform Compucell3D
# CBO CBO_Process CellProcess CellGrowth
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

         # CBO CBO_Process CellProcess CellDeath Necrosis
         # GO:0008219 cell death			   
         #Necrotic Cells
         if cell.type == 3:
            #set growth rate equation
            cell.targetVolume-=0.5
            cell.lambdaSurface=0
	 
      #self.output_file.write("\n")