Autômato Celular para simulação da Esquistossomose em uma região.
Izabela Ralime / Manoel Rodrigues /Ronaldo Ramos
       25 de junho de 2008

((Autômato Celular da em Esquistossomose para Ralime Ramos Rodrigues simulação uma região.Izabela)/(Manoel Ronaldo))

(* 1- ESCOLHA DO DIRETÓRIO*)
(* --Nesta parte há uma interação com o usuário, de onde ele deve dizer que estão os arquivos de dados necessários para o sistema. ---- *)

If[ChoiceDialog["ESTE É O CÓDIGO DE UM AUTÔMATO CELULAR QUE SIMULA A EXPANSÃO DA ESQUISTOSSOMOSE NUMA DETERMINADA ÁREA.

AUTORES:
        IZABELA RALIME
        MANOEL RODRIGUES
        RONALDO RAMOS"],
dados=InputString["ENTRE COM O NOME DO DIRETÓRIO ONDE ESTÁ SALVO O BANCO DE DADOS DO LEVANTAMENTO DE CAMPO:

A pasta deve estar salva no diretório C:"]];
(* ------------------------------------------------------------------------------------- *)

(* 2- MATRIZES AUXILIARES*)

(* Fução que cria matrizes nulas de dimensão (m x n) e (m x n x 3) *)
(* ------------------------------------------------------------------------------------- *)
matrizAuxMolusco[A_]:=(Table[0,{i,Dimensions[A][[1]]},{j,Dimensions[A][[2]]},{k,3}]) 
(* Cria uma matriz de listas de três zeros (cada um indicando um tipo de molusco). As opções Dimensions vão dar as dimensões da matriz de acordo com as dimensões da população A_. O k serve para indicar que devemos preencher cada elemento com 3 zeros. O link dessa função só é gerado, quando a mesma é chamada dentro da função criarMatrizMolusco

OBSERVAÇÃO IMPORTANTE = É UMA MATRIZ TRIDIMENSIONAL, CUJA 3ª DIMENSÃO TEM TAMANHO MÁXIMO IGUAL A 3*)

(* NOVAS MATRIZES A SEREM COMPLETAMENTE IMPLEMENTADAS *)

semEsgoto=Import["C:/"<>dados<>"/esgotoSemTratamento.txt","Table"];

matrizPlotagem[A_] := (Table[0, {i, Dimensions[A][[1]]}, {j, Dimensions[A][[2]]}]);

matrizColecoesAgua[A_]:= (Table[0, {i, Dimensions[A][[1]]}, {j, Dimensions[A][[2]]}]);

matrizEsgoto[A_]:=(Table[0, {i, Dimensions[A][[1]]}, {j, Dimensions[A][[2]]}]);

matrizFinalMolusco[A_]:=(Table[0, {i, Dimensions[A][[1]]}, {j, Dimensions[A][[2]]}]);

(*  matrizRisco[] = ((matrizEsgoto)+(matrizColecoesAgua+)(matrizFinalMolusco))   *)

matrizAux[A_]:=(Table[0,{i,Dimensions[A][[1]]},{j,Dimensions[A][[2]]}]);
(* Cria uma matriz bidimensonal auxiliar de zeros nas dimensões da matriz A_*)

(* ------------------------------------------------------------------------------------- *)

matrizColecoesAgua+matrizEsgoto+matrizFinalMolusco
matrizRisco

(* 3- IMPORTAÇÃO DOS DADOS*)

(* Iportação dos dados a serem utilizados na simulação*)
(* ------------------------------------------------------------------------------------- *)
populacaoSuscetivel=Import["C:/"<>dados<>"/populacaoSucetivel.txt","Table"];
populacaoInfectada=Import["C:/"<>dados<>"/populacaoInfectada.txt","Table"];
populacaoRecuperada=Import["C:/"<>dados<>"/populacaoRecuperada.txt","Table"];
(* As três instruções acima importam e associam para cada variável as matrizes de população suscetível, população infectada e população recuperada *)

criarmatrizMoluscoSadio[A_]:=(a=matrizAuxMolusco[A];For[i=1,i≤Dimensions[A][[1]],i++,a[[i]]=Import["C:/"<>dados<>"/populacaoMoluscoSadio/"<>IntegerString[i,10]<>".txt","Table"]];Return[a])

(* Estra função vai criar a matriz de moluscos de acordo com a seguinte lógica:
1º ==> Teremos associado a variável "a" a matriz de moluscos (já comentada),
  NÃO SEI AO CERTO, VERIFICAR ==> 2º ==> Teremos um loop For que, de acordo com o tamanho da linha da matriz passada vai incrementar a variável "i". Assim em cada linha de "a" será passada a população de molusco sadio??? 
EM resumo.... essa função cria uma nova matriz com a população de moluscos sadios
Vai pegar o tamanho das linhas da matriz A_ e vai criar uma matriz do tipo [tamanho da linha de A][tamanho da linha de A], neste caso é uma matriz do tipo 4x4. Na verdade esta matriz é um vetor de tamanho 4, onde cada elemento é uma matriz do tipo 4x3 (matriz tridimensional). Cada elemento é uma matriz dos arquivos 1.txt, 2.txt, 3.txt e 4.txt que se encontram nos arquivos de populacaoMoluscoSadio *)

populacaoMoluscoSadio=criarmatrizMoluscoSadio[populacaoSuscetivel];
(*Associa a populacaoMoluscoSadio a matriz de moluscos sadios já comentada na linha anterior *)

criarmatrizMoluscoInfectada[A_]:=(a=matrizAuxMolusco[A];For[i=1,i≤Dimensions[A][[1]],i++,a[[i]]=Import["C:/"<>dados<>"/populacaoMoluscoInfectada/"<>IntegerString[i,10]<>".txt","Table"]];Return[a])

(*Faz a mesma coisa que a função criarmatrizMoluscoSadio, o único diferencial é que aqui estamos pegando a população de moluscos infectada dos arquivos de populacaoMoluscoInfectada *)

populacaoMoluscoInfectada=criarmatrizMoluscoInfectada[populacaoSuscetivel];

(*Associa a populacaoMoluscoInfecta a matriz criada de moluscos infectados comentado anteriormente*)

dadosChuva=Import["C:/"<>dados<>"/dadosChuva.txt","Table"];

(*importa e associa a dadosChuva os indices pluviometricos de chuva de aproximadamente 30 em 30 dias, os indices serão associados 2 a 2, tipo:
dia 1 = indice x, dia 15 = indice y, dia 45 = indice z ==> no formato {1, x}, {15, y}, {45, z}. Onde neste caso os valores foram representados por x, y, z, mas os mesmos existem no arquivo dadosChuva.txt*)

areaAlagadaMedia=Import["C:/"<>dados<>"/areaAlagadaMedia.txt","Table"];
(*pega a matriz com as medias das areas alagadas *)

criarmatrizConexaoAgua[A_]:=(a=matrizAuxMolusco[A];For[i=1,i≤Dimensions[A][[1]],i++,a[[i]]=Import["C:/"<>dados<>"/conexao/"<>IntegerString[i,10]<>".txt","Table"]];Return[a])

(* faz uma matriz 4x4 (para este caso),onde cada célula é um vetor de 8 elementos que indica a ligação de agua com os vizinhos 

	ou seja, cria uma matriz tridimensional onde a dimensão z tem oito itens.

Na verdade, esse vetor vai receber em suas linhas o esquema de ligações entre sua vizinhança para cada célula. Percebe-se que existem oito possibilidades para cada célula (vizinhança de Moore), que vão funcionar podendo ter até 5 estados (0, 1, 2, 3 ou 4), onde:

  0 - Não há conexão com a célula
1 - Célula de onde vem o fluxo d'água
2 - Célula para onde vai o fluxo d'água
3 - Fluxo de água vem de uma área externa ao mundo
4 - Fluxo de água vai para uma área externa ao mundo

*)




conexaoAgua=criarmatrizConexaoAgua[populacaoSuscetivel];
(* ------------------------------------------------------------------------------------- *)

(* 4- DRISTIBUIÇÃO DAS CHUVAS *)

(* Dados retirado da estação Goiana *)
(* ------------------------------------------------------------------------------------- *)

interpolacaoDados = Interpolation[dadosChuva];
(* A ordem de interpolação padrão para o mathematica é 3, no entanto este polinomio poderia ser de ordem até 13
Associa os dias com os indices pluviométricos 2 a 2, tipo x associado a y, (x, y) *)

chuvaAreaAlagada[A_,d_]:=(a=matrizAux[A];For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,a[[i,j]]=A[[i,j]]*(interpolacaoDados [d]/200)]];Return[a])
(*a variável d_ são os dias que vão variar no decorrer do ano (de 1 a 365 ) *)
(*no lugar de A_ é utilizado a matriz areaAlagadaMedia, que pega a matriz com as medias das areas alagadas => indices pluviométricos de 30 em 30 dias *)



Plot[interpolacaoDados[x], {x,1,365}, Epilog->Map[Point, dadosChuva],AxesOrigin->{1,0},AxesLabel->{dias,mm³ },ColorFunction->Function[{x,y},Hue[y]]]
(* ------------------------------------------------------------------------------------- *)

(* Grafico ilustrativo da distribuição das chuvas durante o ano *)

(* 5- CRESCIMENTO DOS MOLUSCOS *)

(* Função de Crescimento Populacional:Modelo de Verhulst *)
(* ------------------------------------------------------------------------------------- *)
n[t_,n0_,l_,k_]:=If[(n0+(l-n0) (E^(-k t)))==0,Return[0],(l n0)/(n0+(l-n0) (E^(-k t)))]
(* ------------------------------------------------------------------------------------- *)
(* FUNÇÃO QUE DEFINE CRESCIMENTO POPULACIONAL DOS MOLUSCOS ==> DESEJÁVEL = ESTUDAR OFUNCIONAMENTO 

PROPOSTA ==> SUBSTITUIR A PRIMEIRA CONDIÇÃO "(n0+(l-n0) \!\(\*SuperscriptBox[\(E\), \(\(-k\)\\ t\)]\))==0" POR n0 == 0, DIMINIUNDO O CUSTO COMPUTACIONAL *)

(* função que retorna o número limite de molusco por celula *)
(* ------------------------------------------------------------------------------------- *)
limiteMolusco[i_,j_]:=(areaAlagada[[i,j]]*40);  (* Limite da população em função do ambiente *)
(* SUPOMOS QUE A INTENÇÃO DESTA FUNÇÃO SEJA INDICAR QUE A QUANTIDADE LIMITE DE MOLUSCOS É DE 40 MOLUSCOS POR m² 
  PROVAVELMENTE PARA SER UTILIZADA DENTRO DO MODELO DE VERHULST *)

constanteMolusco=0.02; (* constante de crescimento para ser usada dentro do modelo de VERHLST *)
(* area alagavel de cada celula *)
(* ------------------------------------------------------------------------------------- *)


(* função de crescimento  do molusco *)
(* ------------------------------------------------------------------------------------- *)
(* entrar com a populacao de molusco sadio e infectado *)
crescimentoMolusco[A_]:=(a=matrizAuxMolusco[A];For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,For[k=1,k≤Dimensions[A][[3]],k++,a[[i,j,k]]=n[1,A[[i,j,k]],limiteMolusco[i,j],constanteMolusco]]]];Return[N[a-A]])
(* ------------------------------------------------------------------------------------- *)
(* ESTA FUNÇÃO APLICA A FUNÇÃO DE CRESCIMENTO POPULACIONAL A CADA ELEMENTO DA matrizAuxMolusco RETORNANDO UMA NOVA MATRIZ QUE CONTÉM APENAS O VALOR NUMÉRICO DO CRESCIMENTO PARA CADA ELEMENTO DA MATRIZ QUE FOI PASSADA *)

(* função que retorna a taxa de infecção dos moluscos *)
(* ------------------------------------------------------------------------------------- *)
taxaInfMolusco[i_,j_]:=(If[(300*taxaSaneamento*areaAlagada[[i,j]])==0,0,populacaoInfectada[[i,j]]/(300*taxaSaneamento*areaAlagada[[i,j]])])
(* ------------------------------------------------------------------------------------- *)
(* É POSSÍVEL FAZER UMA MELHOR ORGANIZAÇÃO DOS COMPONENTES DESTE MODELO PARA MELHOR LEGIBILIDADE. OBSERVAR taxaSaneamento*)

(* POR QUE O VALOR 300??? SERIA MAIS UM NÚMERO DE AJUSTE???? *)

(* função de crescimento  do molusco infectado *)
(* ------------------------------------------------------------------------------------- *)
(* entra com populacaoMoluscoSadio e populacaoMoluscoInfectado *)
crescimentoMoluscoInfectado[A_,B_]:=(a=matrizAuxMolusco[A];b=A+B;For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,For[k=1,k≤3,k++,If[B[[i,j,k]]≤b[[i,j,k]]/100,For[w=1,w≤A[[i,j,k]],w++,If[RandomReal[]≤taxaInfMolusco[i,j],a[[i,j,k]]=a[[i,j,k]]+1]]]]]];Return[a])
(* ------------------------------------------------------------------------------------- *)

(*
1 cria uma matriz auxiliar 3d zerada para guardar os novos moluscos infectados.
2 cria laços aninhados para percorrer cada celula da matriz
3 verifica se a populacao infectada de cada celula é inferior a 100%, caso positivo, gera um numero aleatorio e compara com a taxa de infecção do molusco da célula, se menor ou igual adiciona 1 à célula na matriz auxiliar. 
4 vai retornar uma matriz 3d com as quantidades dos novos moluscos infectados de cada espécie. 
*)

(* 6- MOVIMENTAÇÃO DOS MOLUSCOS *)

(* função que simula o deslocamento do molusco pelas celulas *)
(* ------------------------------------------------------------------------------------- *)
taxaMovMoluscoSubir[d_]:=(1/(100*interpolacaoDados[d]))

(* De acordo com o dia vai retornar o valor de 1 dividido por 100 vezes o retorno da função
	segue uma função que vai determinar a taxa para que o molusco suba o rio num determinado dia (período).

 *)


taxaMovMoluscodescer[d_]:=(interpolacaoDados[d]/100000)

(* De acordo com o dia vai retornar o valor da função dividido por 100000 

	segue uma função que vai determinar a taxa para que o molusco desca o rio num determinado dia (período).

A taxa de descer sempre é maior que a de subir.

*)


(* matriz que indica a conexão das celulas com suas vizinhas: 0 - Não possui conexão com a célula, 1->sobe o fluxo, 2->desce o fluxo, 3-> sobe o fluxo p/fora do mundo, 4->desce o fluxo p/fora do mundo *)

escolha1[i_,j_]:=(d=RandomChoice[Position[conexaoAgua[[i,j]],1]][[1]];Which[d==1,{x=i-1,y=j-1},d==2,{x=i-1,y=j},d==3,{x=i-1,y=j+1},d==4,{x=i,y=j-1},d==5,{x=i,y=j+1},d==6,{x=i+1,y=j-1},d==7,{x=i+1,y=j},d==8,{x=i+1,y=j+1}])

(* 
- escolha1[i,j] definirá as coordenadas (x,y) de destino para movimentação do molusco SUBINDO o rio.
- d guardará o número da célula vizinha para onde o molusco se moverá. Vizinhança de Moore:
-- RandomChoice[Position[conexaoAgua[[i,j]],1]][[1]] vai escolher uma das possíveis células vizinhas que tem conexão de água. --

 - as coordenadas de destino (x,y) vão variar de acordo com o valor de d. A referência (i,j) é a célula central. exemplo: d=8 implica x=i+1, y=j+1.
	
*)

(* funções que determina qual deve ser a celula escolhida para descer o fluxo d'água  *)
escolha2[i_,j_]:=(d=RandomChoice[Position[conexaoAgua[[i,j]],2]][[1]];Which[d==1,{x=i-1,y=j-1},d==2,{x=i-1,y=j},d==3,{x=i-1,y=j+1},d==4,{x=i,y=j-1},d==5,{x=i,y=j+1},d==6,{x=i+1,y=j-1},d==7,{x=i+1,y=j},d==8,{x=i+1,y=j+1}])

(* 
- escolha1[i,j] definirá as coordenadas (x,y) de destino para movimentação do molusco DESCENDO o rio.
- d guardará o número da célula vizinha para onde o molusco se moverá. Vizinhança de Moore:
-- RandomChoice[Position[conexaoAgua[[i,j]],1]][[1]] vai escolher uma das possíveis células vizinhas que tem conexão de água. --

 - as coordenadas de destino (x,y) vão variar de acordo com o valor de d. A referência (i,j) é a célula central. exemplo: d=8 implica x=i+1, y=j+1.
	
*)

passeioMolusco[A_,dias_]:=(a=matrizAuxMolusco[A];b=matrizAuxMolusco[A];For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,For[w=1,w≤Dimensions[A][[3]],w++,For[k=1,k≤A[[i,j,w]],k++,
	If[MemberQ[conexaoAgua[[i,j]],1]==True,(* Se tiver conexão de onde vem água *)
		If[RandomReal[]≤taxaMovMoluscoSubir[dias],aa=escolha1[i,j];a[[aa[[1]],aa[[2]],w]]=a[[aa[[1]],aa[[2]],w]]+1;b[[i,j,w]]=b[[i,j,w]]+1, (* Se estiver dentro da taxa de molusco subir de acordo com o dia faz escolha1 *)
			If[MemberQ[conexaoAgua[[i,j]],2]==True, (* Se tiver conexão para onde vai água *)
				If[RandomReal[]≤taxaMovMoluscodescer[dias],aa=escolha2[i,j]; (* Se estiver dentro da taxa do molusco descer o rio de acordo com o dia faz escolha2 *)a[[aa[[1]],aa[[2]],w]]=a[[aa[[1]],aa[[2]],w]]+1;b[[i,j,w]]=b[[i,j,w]]+1], (* ??? adiciona mais um na contagem ??? *)
					If[MemberQ[conexaoAgua[[i,j]],4]==True,(* Se o fluxo de água vai para uma área externa do mundo *)
						If[RandomReal[]≤taxaMovMoluscodescer[dias],b[[i,j,w]]=b[[i,j,w]]+1]]]],
								If[MemberQ[conexaoAgua[[i,j]],3]==True,(* SE o fluxo de água vem de uma área externa do mundo *)
									If[RandomReal[]≤taxaMovMoluscoSubir[dias],b[[i,j,w]]=b[[i,j,w]]+1, (* Se estiver dentro da taxa do molusco subir  naquele dia *)
										If[MemberQ[conexaoAgua[[i,j]],2]==True,(* Se tiver conexão para onde vai a água *)
											If[RandomReal[]≤taxaMovMoluscodescer[dias],aa=escolha2[i,j]; (* Se estiver dentro da taxa do molusco descer *) a[[aa[[1]],aa[[2]],w]]=a[[aa[[1]],aa[[2]],w]]+1;b[[i,j,w]]=b[[i,j,w]]+1]]],
												If[MemberQ[conexaoAgua[[i,j]],2]==True, (* Se tiver conexão para onde vai a água *)
													If[RandomReal[]≤taxaMovMoluscodescer[dias], (* Se estiver dentro da taxa do molusco descer *) aa=escolha2[i,j];a[[aa[[1]],aa[[2]],w]]=a[[aa[[1]],aa[[2]],w]]+1;b[[i,j,w]]=b[[i,j,w]]+1],
													      If[MemberQ[conexaoAgua[[i,j]],4]==True,(* Se o fluxo de água vai para área externa do mundo *)
													       If[RandomReal[]≤taxaMovMoluscodescer[dias],b[[i,j,w]]=b[[i,j,w]]+1]]]]]]]]];Return[{a,b}]) (* Se estiver dentro da taxa do molusco descer *)


(* OBSERVAR ATENTAMETNE O QUE ESTÁ SENDO FEITO ACIMA, PRINCIPALMENTE O PQ ESTES IFs ESTÁO ANINHADOS 
- essa função será responsável pelo movimento dos moluscos, chamada a cada ciclo;
 - recebe a matriz de moluscos infectados e um número de 1 a 365 representando o dia;
 - executa a rotina para cada unidade de molusco infectado, presente em cada célula da matriz 3d de moluscos infectados:

-
*) 
(* ------------------------------------------------------------------------------------- *)

(* 7- TAXAS E CONSTANTES RELATIVAS A POPULAÇÃO HUMANA *)

(* taxas relativas a população humana *)
(* ------------------------------------------------------------------------------------- *)
taxaSaneamento=0.4;
taxaCura=0.005;
taxaMovimentacao=0.8;

costanteInf=208333;
costanteReInf=240000;

(* função que retorna a taxa de infecção da celula *)
(*Fórmula de acordo com o artigo*)
taxaInf[i_,j_]:=If[areaAlagada[[i,j]]*taxaSaneamento* costanteInf==0,0, (((populacaoMoluscoInfectada[[i,j,1]]*400+populacaoMoluscoInfectada[[i,j,2]]*4500+populacaoMoluscoInfectada[[i,j,3]]*2500))/ (areaAlagada[[i,j]]*taxaSaneamento* costanteInf))];

(* Fórumla de acordo com o artigo.*)
(* função que retorna a taxa de reinfecção da celula *)
taxaReInf[i_,j_]:=If[areaAlagada[[i,j]]*taxaSaneamento* costanteInf==0,0,(((populacaoMoluscoInfectada[[i,j,1]]*400+populacaoMoluscoInfectada[[i,j,2]]*4500+populacaoMoluscoInfectada[[i,j,3]]*2500))/(areaAlagada[[i,j]]*taxaSaneamento* costanteReInf))];

(*Média dos tipos de moluscos infectados*)
m= (populacaoMoluscoInfectada[[i,j,1]]*400+populacaoMoluscoInfectada[[i,j,2]]*4500+populacaoMoluscoInfectada[[i,j,3]]*2500)/(areaAlagada[[i,j]]*taxaSaneamento* costanteReInf)
(* glabrsta=4500, straminea=400, tenagophila=2500 *)
(* ------------------------------------------------------------------------------------- *)

(1/(areaAlagada[[5,j]]))0.000010416666666666666` (400 {{{0,0,0},{0,0,0},{0,0,0},{0,0,0}},{{0,0,0},{0,0,0},{0,0,0},{0,0,0}},{{0,0,0},{15,10,0},{0,0,0},{0,0,0}},{{0,0,0},{0,0,0},{0,0,0},{0,0,0}}}[[5,j,1]]+4500 {{{0,0,0},{0,0,0},{0,0,0},{0,0,0}},{{0,0,0},{0,0,0},{0,0,0},{0,0,0}},{{0,0,0},{15,10,0},{0,0,0},{0,0,0}},{{0,0,0},{0,0,0},{0,0,0},{0,0,0}}}[[5,j,2]]+2500 {{{0,0,0},{0,0,0},{0,0,0},{0,0,0}},{{0,0,0},{0,0,0},{0,0,0},{0,0,0}},{{0,0,0},{15,10,0},{0,0,0},{0,0,0}},{{0,0,0},{0,0,0},{0,0,0},{0,0,0}}}[[5,j,3]])

(* 8- INFECÇÃO HUMANA *)



(* função que cria uma matriz com a soma dos humanos infectados apos o passeio *)
(* ------------------------------------------------------------------------------------- *)

(* função que cria uma matriz com a soma dos humanos infectados *)
(* ------------------------------------------------------------------------------------- *)

(*Simplesmente verifica se um valor randomico está dentro da taxa, se estiver adiciona mais um humano infectado para aquela celula*)
somaInfectado[A_]:=(a=matrizAux[A];For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,For[k=1,k≤A[[i,j]],k++,If[RandomReal[]≤taxaInf[i,j],a[[i,j]]=a[[i,j]]+1]]]];Return[a])

(* Vai caminhar na matriz 3D e para cada elemento que esteja dentro da taxa de infectados (taxaInf), vai somar mais um para infectados. *)

(* ------------------------------------------------------------------------------------- *)

(* função que cria uma matriz com a soma dos humanos curados *)
(* ------------------------------------------------------------------------------------- *)
somaCurados[A_]:=(a=matrizAux[A];For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,For[k=1,k≤A[[i,j]],k++,If[RandomReal[]≤taxaCura,a[[i,j]]=a[[i,j]]+1]]]];Return[a])

(* VAi caminhar na matriz 3D e para cada elemento que esteja dentro da taxa de cura (taxaCura), vai somar mais um para os curados *)
(* ------------------------------------------------------------------------------------- *)

(* função que cria uma matriz com a soma dos humanos Reinfectados *)
(* ------------------------------------------------------------------------------------- *)
somaReInfectado[A_]:=(a=matrizAux[A];For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,For[k=1,k≤A[[i,j]],k++,If[RandomReal[]≤taxaReInf,a[[i,j]]=a[[i,j]]+1]]]];Return[a])
(* Vai caminha na matriz 3D e para cada elemento que esteja dentro da taxa de reinfectados (taxaReInf) vai somar mais um para os reinfecatados *)

(* ------------------------------------------------------------------------------------- *)

(* função que provoca a movimentação pela vizinhança (8 vizinhos) *)
(* ==> Vai definir a distância máxima de passeio pela vizinhança de forma aleatória, neste caso a diferença é de no máximo 1 entre vizinhos (nos eixos x e y). *)
(* ------------------------------------------------------------------------------------- *)
passeio[i_,j_,A_]:=(If[i==1,{xmin=1,
xmax=i+1},If[i==Dimensions[A][[1]],{xmin=i-1,xmax=Dimensions[A][[1]]},{xmin=i-1,xmax=i+1}]];If[j==1,{ymin=1,
ymax=j+1},If[j==Dimensions[A][[2]],{ymin=j-1,ymax=Dimensions[A][[2]]},{ymin=j-1,ymax=j+1}]];xx={i,j};yy={i,j};While[ xx==yy,yy={RandomInteger[{xmin,xmax}],RandomInteger[{ymin,ymax}]}];Return[yy])

(* ==> Para cada molusco de uma determinada célula da matriz vai verificar se o mesmo se encontra dentro da taxa de movimentação, caso esteja, vai aplicar a função passeio[i, j, A] e no mesmo bloco vai verificar se está também dentro da taxa de infecção, se estiver, vai incrementar o número de infectados. *)
somaInfectadoPasseio[A_]:=(a=matrizAux[A];For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,For[k=1,k≤A[[i,j]],k++,If[RandomReal[]≤taxaMovimentacao,b=passeio[i,j,A];If[RandomReal[]≤taxaInf[b[[1]],b[[2]]],a[[i,j]]=a[[i,j]]+1]]]]];Return[a])
(* ------------------------------------------------------------------------------------- *)

(* ==> Para cada molusco de uma determinada célula da matriz vai verificar se o mesmo se encontra dentro da taxa de movimentação, caso esteja, vai aplicar a função passeio[i, j, A] e no mesmo bloco vai verificar se está também dentro da taxa de reinfecção, se estiver, vai incrementar o número de reinfectados. *)

(* função que cria uma matriz com a soma dos humanos Reinfectados apos o passeio *)
(* ------------------------------------------------------------------------------------- *)
somaReInfectadoPasseio[A_]:=(a=matrizAux[A];For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,For[k=1,k≤A[[i,j]],k++,If[RandomReal[]≤taxaMovimentacao,b=passeio[i,j,A];If[RandomReal[]≤taxaReInf[b[[1]],b[[2]]],a[[i,j]]=a[[i,j]]+1]]]]];Return[a])
(* ------------------------------------------------------------------------------------- *)

(* 9- NOVAS GERAÇÕES *)


(* função que calcula a razão infectatos por população total *)

(* ------------------------------------------------------------------------------------- *)
porcentagemPopulacao[A_,B_,C_]:=(a=matrizAux[A];For[i=1,i≤Dimensions[A][[1]],i++,For[j=1,j≤Dimensions[A][[2]],j++,If[A[[i,j]]+B[[i,j]]+C[[i,j]]==0,a[[i,j]]=0,a[[i,j]]=A[[i,j]]/(A[[i,j]]+B[[i,j]]+C[[i,j]])]]];Return[IntegerPart[10*N[a]]])
(* ------------------------------------------------------------------------------------- *)

(* função que calcula a proxima geração *)
(* ------------------------------------------------------------------------------------- *)
anos=0;

(* Função principal que será controlada no Manipulate, vai controlar todo o ciclo de gerações *)
(* No primeiro caso (dia 1) e no último caso (dia 365) retorna a porcentagem de população (infectada, suscetível e recuperada) em caso diferente vai começar a utilizar todas as funções já estudadas, utilizando sempre o dia atual (de onde será retirado o indice pluviométricos e tudo de acordo com as outras funções*)


(*ESTUDO*)

matrizTotalInfectados[A_]:=
(a=matrizPlotagem[A];
	For[i=1,i≤Dimensions[A][[1]],i++,
		For[j=1,j≤Dimensions[A][[2]],j++,
			For[k=1,k≤3,k++,
				a[[i, j]]= a[[i,j]] + A[[i, j,k]]
			]
		]
	];
Return[a])
(* ------------------------------------------------------------------------------------- *)

totalInfectados=matrizTotalInfectados[populacaoMoluscoInfectada]

semEsgoto=Import["C:/"<>dados<>"/esgotoSemTratamento.txt","Table"]

localAreaAlagada[A_]:=(a=matrizPlotagem[A];
						For[i=1, i≤Dimensions[A][[1]], i++,
							For[j=1, j≤Dimensions[A][[2]], j++,
								If[A[[i, j]] == 0,a[[i, j]] = 0, a[[i, j]] = 1]]];
						Return[a])


flagAreaAlagada = localAreaAlagada[areaAlagada]

level1 = 20;
level2 = 45;
level3 = 70;


NextGen = {{_,x_/;0≤ x< level1}->0,
		{_,x_/;level1≤ x< level2}->1,
		{_,x_/;level2≤ x< level3}->2,
		{_,x_/;level3≤ x}->3};


(* convertMatriz[A_]:=(a=matrizPlotagem[A];
					For[i=1, i≤Dimensions[A][[1]], i++,
						For[j=1, j≤Dimensions[A][[2]], j++,
							]];
					Return[a])  *)

EightNeighbors [earth_]:=earth+(RotateRight /@ earth) + (RotateLeft /@ earth) + RotateRight /@(RotateRight[earth])+ (RotateRight[earth])+ RotateLeft /@(RotateRight[earth])+ RotateRight /@(RotateLeft[earth])+ (RotateLeft[earth])+ RotateLeft /@(RotateLeft[earth]);

EightNeighbors[populacaoInfectada]


risco[A_] := (semEsgoto+totalInfectados+flagAreaAlagada+EightNeighbors[A])

CarneDeVaca[a_] := (Transpose/@ Transpose[{ a,risco[a]}])/. NextGen;

RodaLaVaca[matriz_, n_]:= If[n == 0, matriz, RodaLaVaca[CarneDeVaca[matriz],n-1]];



proximaGeracao[a_] := (If[a==1,Return[RodaLaVaca[populacaoInfectada, a]],areaAlagada=chuvaAreaAlagada[areaAlagadaMedia,a];populacaoMoluscoSadio=populacaoMoluscoSadio+crescimentoMolusco[populacaoMoluscoSadio]+crescimentoMolusco[populacaoMoluscoInfectada];moluscoInfectadoFuncao=crescimentoMoluscoInfectado[populacaoMoluscoSadio,populacaoMoluscoInfectada];populacaoMoluscoInfectada=populacaoMoluscoInfectada+moluscoInfectadoFuncao;populacaoMoluscoSadio=populacaoMoluscoSadio-moluscoInfectadoFuncao;moluscoSadioPasseio=passeioMolusco[populacaoMoluscoSadio,a];moluscoinfectadoPasseio=passeioMolusco[populacaoMoluscoInfectada,a];populacaoMoluscoSadio=populacaoMoluscoSadio+moluscoSadioPasseio[[1]]-moluscoSadioPasseio[[2]];populacaoMoluscoInfectada=populacaoMoluscoInfectada+moluscoinfectadoPasseio[[1]]-moluscoinfectadoPasseio[[2]];populacaoInfPasseio=somaInfectadoPasseio[populacaoSuscetivel];populacaoReInfPasseio=somaReInfectadoPasseio[populacaoRecuperada];populacaoSuscetivel=populacaoSuscetivel-populacaoInfPasseio;populacaoRecuperada=populacaoRecuperada-populacaoReInfPasseio;populacaoInfectada=populacaoInfectada+populacaoInfPasseio+populacaoReInfPasseio;infectados=somaInfectado[populacaoSuscetivel];reinfectados=somaReInfectado[populacaoRecuperada];
curados=somaCurados[populacaoInfectada];populacaoInfectada=populacaoInfectada+infectados+reinfectados-curados;populacaoSuscetivel=populacaoSuscetivel-infectados;populacaoRecuperada=populacaoRecuperada+curados-reinfectados;If[a==365,anos++];Return[RodaLaVaca[populacaoInfectada, a]]])

(* 10- RELATÓRIO *)
(* função que gera um relatório por celula *)
relatorio[i_,j_]:=(Print[StringForm["população sucetível =          ``",populacaoSuscetivel[[i,j]]]]Print[StringForm["população recuperada =         ``",populacaoRecuperada[[i,j]]]]Print[StringForm["população infectada =          ``",populacaoInfectada[[i,j]]]]Print[StringForm["porcentagem infectada =        ``%",NumberForm[100*N[(populacaoInfectada[[i,j]]/(populacaoSuscetivel[[i,j]]+populacaoRecuperada[[i,j]]+populacaoInfectada[[i,j]]))],{2,2}]]]Print[StringForm["população molusco sadio =      ``",IntegerPart[populacaoMoluscoSadio][[i,j]]]]Print[StringForm["população molusco infectada =  ``",populacaoMoluscoInfectada[[i,j]]]]Print[StringForm["área alagada =                 ``m²",NumberForm[areaAlagada[[i,j]],{2,2}]]]Print[StringForm["taxa de infecção =             ``%",NumberForm[100*taxaInf[i,j],{2,2}]]]Print[StringForm["taxa de Reinfecção =           ``%",NumberForm[100*taxaReInf[i,j],{2,2}]]]Print[StringForm["taxa infc. molusco =           ``%",NumberForm[100*taxaInfMolusco[i,j],{2,2}]]];Return["Fim do relatório"])

(* 11- INTERFACE GRÁFICA *)

MostrarAreaInfectada[a_]:=Show[ArrayPlot[a,Mesh->True,ColorRules->{0-> LightYellow,1-> Yellow,2->LightRed,3->Red}]]

Import["C:/"<>dados<>"/Drawing1.bmp"]

(* função que mostra as gerações uma a uma *)
(* ------------------------------------------------------------------------------------- *)
Manipulate[MostrarAreaInfectada[proximaGeracao[dias]],{dias,1,365,1},"diminuir o slower"]

relatorio[1,5]

populacaoSuscetivel

populacaoInfectada

populacaoRecuperada

IntegerPart[populacaoMoluscoSadio]

populacaoMoluscoInfectada

anos

2
TableForm[conexaoAgua[]]
Position[conexaoAgua[[3,3]], 1]

RandomChoice[Position[conexaoAgua[[3,3]],1]][[1]]

escolha1[3,3]

m

m= (populacaoMoluscoInfectada[[i,j,1]]*400+populacaoMoluscoInfectada[[i,j,2]]*4500+populacaoMoluscoInfectada[[i,j,3]]*2500)/(areaAlagada[[i,j]]*taxaSaneamento* costanteReInf)

populacaoMoluscoInfectada

(*  EM ESTUDOOOOOOOO *)

areaAlagada

semEsgoto

{{0,1,0,1},{0,1,0,1},{0,1,1,1},{0,0,1,0}}

{{28,38,22,27},{71,81,72,18},{81,125,84,34},{64,79,72,24}}

{{28,38,22,27},{71,81,72,18},{81,125,84,34},{64,79,72,24}}
CarneDeVaca[populacaoInfectada]

RodaLaVaca[populacaoInfectada, 2]

RodaLaVaca[populacaoInfectada, 3]

ShowInfectedArea[a_]:=Show[ArrayPlot[a,Mesh->True,ColorRules->{0-> White,1->LightGray,2->Gray,3->Black}]]


Manipulate[ShowInfectedArea[RodaLaVaca[infectedarea,x]],{x,0,20, 1}]