source("fdim.R")
source("ca1.R")

#fractal dimensions of lattices generated by cellular automata


dim= 50
iffig=1
nostep=20
animation=0
type=1

#species 1
imi1=1
ext1=0.7
long1=0.0
ini1=0.1
druh1 = c(imi1,ext1,long1,ini1)

#one species
nospec= 1
specpara=rbind(NULL,druh1)
cc=as.matrix(1,ncol=1)

a=ca(nospec, specpara, cc,dim,type, 1, nostep, animation,1)
fdim(a$map,1)

species 2
imi2=0.7
ext2=0.05
long2=0.0
ini2=0.1
druh2 = c(imi2,ext2,long2,ini2)

#two species
nospec=2
specpara=rbind(druh1,druh2)
cc = matrix(c(1,1,0,1),nrow=nospec)

a=ca(nospec, specpara, cc,dim,type, 1, nostep, animation,1)

spec=1 #select species for getting its fractal dimension
specmap = matrix(as.numeric(a$map == spec),nrow=nrow(a$map),ncol=ncol(a$map))
fdim(specmap,1)


#finding fractal structure of the species at the verge of extinction

ext1=0.4
long1=0.0
ini1=0.1
nospec= 1
cc=as.matrix(1,ncol=1)

fd = NULL
nl = NULL
rge = seq(0.08,0.2,by=0.01)
for (imi1 in rge) {
  druh1 = c(imi1,ext1,long1,ini1)
  specpara=rbind(NULL,druh1)
  a=ca(nospec, specpara, cc,dim,type, 0, nostep, animation,1)
  xx=fdim(a$map,1)
  fd = c(fd,xx[1])
  nl = c(nl,xx[2])
}

plot(rge,fd,ylab="Boxcounting dimension",y="Natality per cell",pch=19)
plot(rge,nl,ylab="Quadratic term",y="Natality per cell",pch=19)
matplot(rge,cbind(fd,nl2),pch=c(19,21),type="b")