function [med, std]= brownie ()

// Movimento browniano:
// dizemos que a trajetoria da variavel aleatoria X
// e' um movimento browniano quando:
// dado x(0),
// x(i+1)=x(i)*(1+u) com probabilidade pu;
// x(i+1)=x(i)*(1-d) com probabilidade pd;
// x(i+1)=x(i) com probabilidade 1-pu-pd.
// u, d > 0.

// No programa vamos decidir se X sobe, desce ou fica na
// mesma baseados em uma variavel aleatoria r com distri-
// buicao uniforme no intervalo [0,1]:

//              pu               pd              1-pu-pd
//        |<------------>|<------------------>|<--------->|
//r:      0--------------|--------------------|-----------1
//        |     pu       |                    |
//        |------------->|                    |
//        |              pu+pd                |
//        |---------------------------------->|

//                  se 0<=r<=pu entao X sobe
//                  se pu<r<=(pu+pd) entao X desce
//                  caso contrario, X fica na mesma.

// Esta versao do programa simulara' a compra de uma opcao.
// Os parametros a serem fornecidos pelo usuario serao:
//     u, d, pu, pd = ja descritos acima
//     N = dia do exercicio da opcao
//     K = preco de exercicio
//     B = limite de knock-out
//     M = numero de simulacoes (trajetorias)

confirma = " "
while confirma<>"S" do
  parametros = input("Digite u, d, pu, pd, N, K, B, M separados por virgulas");

  u = parametros(1);
  d = parametros(2);
  pu = parametros(3);
  pd = parametros(4);
  N = int(parametros(5));
  K = parametros(6);
  B = parametros(7);
  M = parametros(8);

  printf("u = %f", u);
  printf("d = %f", d);
  printf("pu = %f", pu);
  printf("pd = %f", pd);
  printf("N = %d", N);
  printf("K = %f", K);
  printf("B = %f", B);
  printf("M = %d", M);

  confirma = input("Confirma (S/N)?","string");
  confirma = strsubst(confirma,"s","S");
end

rand('Uniform', timer());
				// Abre arquivo de resultados parciais
arqsai=file("open", "exemp3.out", "unknown");

med = 0.0;
var = 0.0;
for j=1:M                        // trajetorias
   x = zeros(N);
   x(1) = 1.0;
   i = 2;
   while i <= N & x(i-1) < B do  // passos
      r=rand(1,1);
      if (r<pu) then
         x(i)=x(i-1)*(1+u);
      elseif (r<(pu+pd)) then
         x(i)=x(i-1)*(1-d);
      else
         x(i)=x(i-1);
      end
      i = i + 1;
   end
   i = i-1

   if i=N & x(i) > K then
      g = x(i)-K;
   else
      g = 0.0;
   end
   fprintf(arqsai,"i = %d    x(i) = %f    ganho = %f", i, x(i), g);

   med = med + g;
   var = var + g^2;
end
med=med/M;
var=var/M - med^2;
std = sqrt(var);

printf("Media = %f", med);
printf("Desvio Padrao = %f", std);

file("close",arqsai);	// Fecha arquivo
