program FRACTAL;
{$E+}{$N+}
uses crt,grhacks,h3,r3,r4,r4x4,graph;

const  xmin=-100;
       xmax=100;
       ymin=-100;
       ymax=100;
       maxbrilho=255;

type  ponto_do_fractal=record
                          t:real;
                          c:h3.point
                        end;
      ponto_de_luz=record
                      coord:h3.point;
                     brilho:r3.T
                   end;

var n,L,v,contx,conty:word;
    lin,col:integer;
    x,y:real;
    cor:r3.T;
    anota:array[1..500] of ponto_do_fractal;
    PL:array[1..10] of ponto_de_luz;
    M:array[1..3] of r4x4.T;
    f,o,u,CE,PE,c,pef:h3.point;
    mtp:h3.pmap;    {matriz de transf. de perspectiva}
    ch,nome:string;
    arq:file of char;

(*******************************************************************)

function quad(x:double): double;
begin
  quad:=x*x
end;

(*******************************************************************)

function min(x,y:double): double;
begin
  if x<y then min:=x else min:=y
end;

(*******************************************************************)

procedure obtem_coord_cartesianas(p:h3.point);

var i:byte;

begin
  for i:=0 to 3 do p.c[i]:=p.c[i]/p.c[0]
end;

(*******************************************************************)

procedure cria_arquivo;

var c:char;

begin
  assign(arq, nome);
  rewrite(arq);
  c:='P'; write(arq, c);
  c:='6'; write(arq, c);
  c:=' '; write(arq, c);
  c:='2'; write(arq, c);
  c:='7'; write(arq, c);
  c:='6'; write(arq, c);
  c:=' '; write(arq, c);
  c:='2'; write(arq, c);
  c:='7'; write(arq, c);
  c:='6'; write(arq, c);
  c:=' '; write(arq, c);
  c:='2'; write(arq, c);
  c:='5'; write(arq, c);
  c:='5'; write(arq, c);
  c:=' '; write(arq, c);
end;

(*******************************************************************)

procedure grava_cor_do_pixel(cor:r3.T);

var r,g,b:byte; c:char;

begin
  r:=round(maxbrilho*cor[0]); c:=chr(r); write(arq,c);
  g:=round(maxbrilho*cor[1]); c:=chr(g); write(arq,c);
  b:=round(maxbrilho*cor[2]); c:=chr(b); write(arq,c);
end;

(*******************************************************************)

procedure le_os_dados;

var x1, y1, z1, x2, y2, z2, x3, y3, z3:double;
    st:string;
    arq:text;
    i,k:word;

begin
  h3.mk_point(1.0,  0.0, 0.0, 0.0, CE);
  h3.mk_point(1.0, xmax, 0.0, 0.0, PE);
  h3.mk_point(1.0,  0.0, 0.0, 0.0, c);
  h3.mk_point(1.0, 10.0, 0.0, 0.0, pef);

  write('Qual o nome do arquivo : '); readln(st);
  assign(arq, st);
  reset(arq);

  readln(arq, n);
  for k:=1 to n
  do begin
       readln(arq);
       for i:=0 to 3 do readln(arq, M[k,i,0], M[k,i,1], M[k,i,2], M[k,i,3]);
     end;

  readln(arq); readln(arq, L); readln(arq);
  for k:=1 to L do
    readln(arq, PL[k].coord.c[0], PL[k].coord.c[1],
                PL[k].coord.c[2], PL[k].coord.c[3],
                PL[k].brilho[0],  PL[k].brilho[1], PL[k].brilho[2]);

  close(arq);

  writeln ; writeln('Entre com obs, foco e u.');
  readln(x1, y1, z1, x2, y2, z2, x3, y3, z3);
  h3.mk_point(1.0, x1, y1, z1, o);
  h3.mk_point(1.0, x2, y2, z2, f);
  h3.mk_point(1.0, x3, y3, z3, u);

  writeln ; write('Entre com nome do arquivo que sera gerado: ');
  readln(nome);

end;

(*******************************************************************)

function muito_proximos(p,q:h3.point) : boolean;

begin
  if (abs(p.c[1] - q.c[1]) < 1.0) and
     (abs(p.c[2] - q.c[2]) < 1.0) and
     (abs(p.c[3] - q.c[3]) < 1.0)
  then muito_proximos:=true
  else muito_proximos:=false
end;

(*******************************************************************)

procedure anota_ponto(x:real; c:h3.point);
begin
  v:=v+1; anota[v].t:=x;
  h3.mk_point(c.c[0], c.c[1], c.c[2], c.c[3], anota[v].c)
end;

(*******************************************************************)

function minimo_da_lista : word;

var  i,m:word; d,d1:double;

begin
  if v=0
  then minimo_da_lista:=0
  else begin
         m:=1;
         for i:=2 to v do if anota[i].t < anota[m].t then m:=i;
         minimo_da_lista:=m
       end
end;

(*******************************************************************)
{
procedure vetor_unitario(a,b:h3.point; var res:h3.point);

var c:double; i:word;

begin
  res.c[0]:=1.0;
  for i:=1 to 3 do res.c[i]:=b.c[i] - a.c[i];
  c:=sqrt(1.0/(quad(res.c[1]) + quad(res.c[2]) + quad(res.c[3])));
  for i:=1 to 3 do res.c[i]:=c*res.c[i];
end;
}
(*******************************************************************)

procedure obtem_ponto_eq_parametrica(p1, p2:h3.point; t:real;
                                    var res:h3.point);
var i:word;

begin
  res.c[0]:=1.0;
  for i:=1 to 3 do res.c[i]:=p1.c[i] + t*(p2.c[i] - p1.c[i])
end;

(*******************************************************************)

procedure resolve_eq_2_gr(a, b, c:double; var nsol:byte; var x1, x2:double);

var delta:double;

begin
  delta:=b*b - 4*a*c;
  if delta < 0
  then nsol:=0
  else if delta = 0
       then begin
              nsol:=1 ; x1 := (-b)/(2*a)
            end
       else begin
              nsol:=2 ;
              x1 := ((-b) + sqrt(delta))/(2*a) ;
              x2 := ((-b) - sqrt(delta))/(2*a)
            end
end;

(*******************************************************************)

function intercepta( c, pr, o, p:h3.point; var t:real) : boolean ;

var r,a,b,d,x1,x2:double;  nsol:byte;

begin
  r:=h3.dst(pr, c);

  a := quad(p.c[1]-o.c[1]) + quad(p.c[2]-o.c[2]) + quad(p.c[3]-o.c[3]) ;
  b := 2*((p.c[1]-o.c[1])*(o.c[1]-c.c[1]) +
          (p.c[2]-o.c[2])*(o.c[2]-c.c[2]) +
          (p.c[3]-o.c[3])*(o.c[3]-c.c[3]));
  d := quad(o.c[1]-c.c[1]) + quad(o.c[2]-c.c[2]) + quad(o.c[3]-c.c[3]) - r*r;

  resolve_eq_2_gr(a, b, d, nsol, x1, x2);

  if (nsol=0) or ((nsol=1) and (x1<=0)) or ((nsol=2) and (x1<=0) and (x2<=0))
  then intercepta:=false
  else begin
         intercepta:=true;
         if nsol=1
         then  t:=x1
         else if x1 <= 0
              then t:=x2
              else if x2<=0
                   then t:=x1
                   else t:=min(x1,x2)
       end
end;

(*******************************************************************)

procedure varre_o_fractal(CE, PE, o, p, c, pef:h3.point);

var CE1, PE1, c1, pef1:h3.point; i:word;
    t:real;

begin
{  if intercepta(CE, PE, o, p, t)
  then}
       begin
         if intercepta(c, pef, o, p, t) then anota_ponto(t,c);
         if not muito_proximos(pef, c)
         then for i:=1 to n do
              begin
                r4x4.map_row( CE.c, M[i],  CE1.c);
                r4x4.map_row( PE.c, M[i],  PE1.c);
                r4x4.map_row(  c.c, M[i],   c1.c);
                r4x4.map_row(pef.c, M[i], pef1.c);
                obtem_coord_cartesianas(CE1);
                obtem_coord_cartesianas(PE);
                obtem_coord_cartesianas(c1);
                obtem_coord_cartesianas(pef1);
                varre_o_fractal(CE1, PE1, o, p, c1, pef1)
              end
       end
end;

(*******************************************************************)

function enxerga(org, dest:h3.point): boolean;

var q:word; p1:h3.point;

begin
  v:=0;
  varre_o_fractal(CE, PE, org, dest, c, pef);
  q:=minimo_da_lista;
  if q=0 then enxerga:=true
         else begin
                obtem_ponto_eq_parametrica(org, dest, anota[q].t, p1);
                if muito_proximos(p1, dest) then enxerga:=true
                                            else enxerga:=false
              end
end;

(*******************************************************************)

procedure calcula_cor_do_ponto(q, centro:h3.point; var res:r3.t);

var cor_prov,cori,f,n:r3.T; fator:double; i:word;
    f1,n1:h3.point;

begin
  cor_prov[0]:=0.0; cor_prov[1]:=0.0; cor_prov[2]:=0.0;
  for i:=1 to L do
   if enxerga(PL[i].coord, q)
   then begin
          h3.dir(q, PL[i].coord, f);
          h3.dir(centro, q, n);
          fator:=r3.dot(f, n);
          if fator > 0.0 then begin
                                r3.scale(PL[i].brilho, fator, cori);
                                r3.add(cor_prov, cori, cor_prov)
                              end

       end;
  res[0]:=cor_prov[0]; res[1]:=cor_prov[1]; res[2]:=cor_prov[2]
end;

(*******************************************************************)

procedure calcula_cor_do_pixel(x, y:real; var res:r3.T);

var p,res1,p1:h3.point; q:word; t:real;

begin
  h3.mk_point(1.0, x, y, 0.0, p);
  r4x4.map_row(p.c, mtp.inv, p.c);
  obtem_coord_cartesianas(p);
  v:=0;
  varre_o_fractal(CE, PE, o, p, c, pef);
  q:=minimo_da_lista;
  if q=0 then begin
                res[0]:=0.0; res[1]:=0.0; res[2]:=0.0
              end
         else begin
                obtem_ponto_eq_parametrica(o, p, anota[q].t, p1);
                calcula_cor_do_ponto(p1, anota[q].c, res)
              end

{
  if intercepta(CE, PE, o, p, t)
  then begin
         v:=0; anota_ponto(t, CE);
         obtem_ponto_eq_parametrica(o, p, anota[1].t, p1);
         calcula_cor_do_ponto(p1, anota[1].c, res)
       end
  else begin
         res[0]:=0.0; res[1]:=0.0; res[2]:=0.0
       end
}
end;

(*******************************************************************)

procedure plota_o_pixel(x, y:real; var cor:r3.T);

var i:word;  c:real;

begin
  for i:=0 to 2 do if cor[i]<0.0 then cor[i]:=0.0;
  for i:=0 to 2 do if cor[i]>1.0 then cor[i]:=1.0;
  c:=0.0;
  for i:=0 to 2 do c:=c + cor[i];
  c:=c/3;
  grhacks.gr_set_color(c);
  grhacks.gr_draw_point(x,y);
end;

(*******************************************************************)
(*******************************************************************)

begin
  le_os_dados;
  h3.persp_map(o, f, u, mtp);

  cria_arquivo;

  grhacks.gr_begin_drawing(xmin, xmax, ymin, ymax);
  grhacks.gr_set_color(0.5);
  grhacks.gr_paint_plane;

  y:=ymax; lin:=grhacks.gr_y(y); conty:=0;
  while y >= -15
  do begin

{       if conty=0 then begin}

       x:=-15; col:=grhacks.gr_x(x); contx:=0;
       while x <= xmax
       do begin
{            contx:=contx+1;}
            if (y <  17*x/18 + 15) and (y > (18/17)*(x - 15))
            then calcula_cor_do_pixel(x,y,cor)
            else begin
                   cor[0]:=0.0; cor[1]:=0.0; cor[2]:=0.0
                 end;

            plota_o_pixel(x,y,cor);

            grava_cor_do_pixel(cor);

            x:=x+0.4;
            while (x<=xmax) and (col=grhacks.gr_x(x)) do x:=x+0.4;
            col:=grhacks.gr_x(x)
          end;
{                                 end;
      conty:=conty + 1;               }

       y:=y-0.4;
       while (y>=-15) and (lin=grhacks.gr_y(y)) do y:=y-0.4;
       lin:=grhacks.gr_y(y)
     end;

 closegraph;
 close(arq);
{ writeln('linhas:',conty,'  colunas:',contx);}

  read(ch);
end.
