#include <iostream.h>

#include "TRandom.h"
#include "TPad.h"
#include "TCanvas.h"
#include "TProfile.h"
#include "TFrame.h"

int teorico1_new()
{

   //Inizializzazione//
   gROOT->Reset();
   TCanvas* c1 = new TCanvas("c1","Modello preda/predatore",200,10,800,800);
   c1->SetFillColor(18);
   TPad* pad1 = new TPad("pad1","Prede",0.05,0.50,0.95,0.95,21);
   TPad* pad2 = new TPad("pad2","Predatori",0.05,0.05,0.95,0.45,21);
   pad1->Draw();
   pad2->Draw();
   pad1->cd();

   //TFile *hfile = gROOT->FindObject("output.root"); if (hfile) hfile->Close();
   //hfile = new TFile("output.root","RECREATE","Conviene averlo per analisi off-line");
 
   const int maxtime = 100;
   /*   const float a = 0.50; 
   const float b = 0.60; 
   const float c = 0.7143;*/
   const float am = 0.01; 
   //   Double_t par[4]={0.7, 0.007 , 0.5 , 0.1 };
   //   Double_t par[4]={0.7, 0.00007 , 0.5 , 0.1 };
   Double_t par[4]={0.7, 0.00003 , 0.5 , 1 };
   int T = maxtime/am ;

   float nprede = 21000;
   float npredatori = 49000;

   float nprede_new=nprede;
   float npredatori_new=npredatori;

   float max1=nprede;
   float max2=npredatori;
   
   TProfile *hprede = new TProfile("hprede","N. prede in funzione del tempo",T,0,T,0,nprede_new*4);
   TProfile *hpredatori = 
   new TProfile("hpredatori","N. predatori in funzione del tempo",T,0,T,0,npredatori_new*4);

   gRandom->SetSeed();
 
   hprede->SetXTitle("cicli");
   hprede->SetYTitle("prede");
   hpredatori->SetXTitle("cicli");
   hpredatori->SetYTitle("predatori");

   hprede->SetMarkerStyle(2);
   hpredatori->SetMarkerStyle(2);
   
   pad1->cd();
   pad1->SetGridx();
   pad1->SetGridy();
   pad1->GetFrame()->SetBorderMode(-1);
   pad1->GetFrame()->SetBorderSize(5);

   pad2->cd();
   pad2->SetGridx();
   pad2->SetGridy();
   pad2->GetFrame()->SetBorderMode(-1);
   pad2->GetFrame()->SetBorderSize(5);
   
   //Fine della inizializzazione//

   //Modello standard//

   for ( Int_t t=1; t<T; t=t+1) 
   {

      // Controlla se il massimo e' stato raggiunto, e in caso aggiornalo //
      
      if (nprede>max1)  
      {
         max1=nprede;
      }
      if (npredatori>max2) 
      {
         max2=npredatori;
      }


      //Riempi gli istogrammi//
    
      hprede->Fill(t,nprede);
      hpredatori->Fill(t,npredatori);
      cout
        << "prede: " << nprede
        << ", predatori: " << npredatori
        << ", time: " << t << " " << maxtime
        << endl;

      //Visualizzazione su grafico punto per punto//

      /*      pad1->cd();
      hprede->Draw("hist p");
      pad1->Update();

      pad2->cd();
      hpredatori->Draw("hist p");
      pad2->Update(); */
      
      //Prima di passare al prossimo ciclo aggiorno i valori delle variabili//

      /*      nprede_new=nprede+am*(nprede-a*nprede*npredatori);
      npredatori_new=npredatori+am*(b*nprede*npredatori-c*npredatori);*/
      nprede_new=nprede+am*(par[0]*nprede-par[1]*nprede*npredatori);
      npredatori_new=npredatori+am*(par[3]*par[1]*nprede*npredatori-par[2]*npredatori);
      nprede=nprede_new;
      npredatori=npredatori_new;

   }
   pad1->cd();
   hprede->SetMaximum(1.25*max1);
   hprede->SetMinimum(-1);
   hprede->Draw("hist p");
   //   pad1->Update();
   
   pad2->cd();
   hpredatori->SetMaximum(1.25*max2);
   hpredatori->SetMinimum(-1);
   hpredatori->Draw("hist p");
   //   pad2->Update(); 
   return 0; 
 }