/* vide {rdo_sampling.h} */
#define rdo_sampling_C_COPYRIGHT "Copyright © 2008 Danillo Pereira e J. Stolfi, UNICAMP"
/* Last edited on 2008-07-22 16:25:50 by stolfi */ 

#define _GNU_SOURCE
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <assert.h>

#include <frgb.h>
#include <argparser.h>
#include <jsrandom.h>

#include "rdo_geometry.h"
#include "rdo_solid.h"
#include "rdo_finish.h"
#include "rdo_site.h"
#include "rdo_lighting.h"
#include "rdo_camera.h"
#include "rdo_radiosity.h"
#include "rdo_sampling.h"

void GeraSitiosPorProjecao
  ( ListaDeObjetos *obs, 
    Camera *cam, 
    int nh,
    int nv,
    ListaDeSitios *sts
  )
  {
    Ponto3d ctr;     //centro da imagem no espaço.
    Vetor3d r, s, t; //vetores unitários dos eixos da câmera no espaço.
    DeterminaEixosDoFilme(cam, &ctr, &r, &s, &t);
    
    //calcula espaçamento médio entre os raios no plano da imagem:
    double nr = hypot((double)nh, (double)nv);
    double passo = cam->raio/nr;

    //percorre a grade de pontos no plano da imagem:
    int nold = sts->ne;  //número original de sítios em {sts}
    int ntot = ntot;     //numero total de sítios em {sts}.
    int kh, kv;
    for(kh = 0; kh < nh; kh++)
      for(kv = 0; kv < nv; kv++)
        {
          //escolhe a posição {(h,v)} do raio no sistema de coordenadas da câmera:
          double h = (kh - ((double)nh)/2 + 0.1 + 0.8*drandom())*passo;
          double v = (((double)nv)/2 - kv + 0.1 + 0.8*drandom())*passo;
          
          //calcula a intersecção {p} do raio com o plano da imagem, no espaço:
          Ponto3d p = ctr;
          r3_mix_in(+h, &r, &p);
          r3_mix_in(-v, &s, &p);

          //determina a direcao do raio:
          Vetor3d dir = DirecaoDePontoParaPonto(&(cam->olho), &p);

          //traça o raio e armazena o resultado:
          Sitio sR;
          double dR;
          SitioDeInterseccaoRaioObjetos(obs, &(cam->olho), &dir, &sR, &dR);
          if (isfinite(dR) && (sR.iob >= 0))
            { ListaDeSitios_expand(sts, ntot);
              sts->e[ntot] = sR;
              ntot++; 
            } 
        }
    fprintf(stderr, "gerados %d sítios por projecao (total %d)\n ", ntot - nold, ntot);
    ListaDeSitios_trim(sts, ntot);
  }

void GeraSitiosSecundarios(ListaDeObjetos *obs, ListaDeSitios *sts, int num_raios)
  { int i; 
    //tenta gerar {num_raios} filhos de cada sítio existente:
    int nold = sts->ne;  //número original de sítios em {sts}
    int ntot = ntot;     //numero total de sítios em {sts}.
    for(i = 0; i < nold; i++) 
      { int k;
        for (k = 0; k < num_raios; k++)
          { Sitio sR = GeraSitioSecundario(obs, &(sts->e[i]));
            if (sR.iob >= 0)
              { ListaDeSitios_expand(sts, ntot);
                sts->e[ntot] = sR;
                ntot++; 
              }
          }
      }
    fprintf(stderr, "gerados %d sítios derivados (total %d)\n ", ntot - nold, ntot);
    ListaDeSitios_trim(sts, ntot);
  }

ListaDeSitios GeraSitiosDosPixeis
  ( ListaDeObjetos *obs,
    Camera *cam, 
    Filme *filme
  )
  {
    int nh = filme->sz[1];
    int nv = filme->sz[2];
    ListaDeSitios sts = ListaDeSitios_new(nh*nv);

    Ponto3d ctr;     //centro da imagem no espaço.
    Vetor3d r, s, t; //vetores unitários dos eixos da câmera no espaço.
    DeterminaEixosDoFilme(cam, &ctr, &r, &s, &t);
    
    //calcula espaçamento médio entre os raios no plano da imagem:
    double nr = hypot((double)nh, (double)nv);
    double passo = cam->raio/nr;

    //percorre os pixeis da imagem:
    int ntot = 0;  //contador do numero de sítios de pixeis:
    int ih, iv;    // Indices do pixel
    for(ih = 0; ih < nh; ih++)
      for(iv = 0; iv < nv; iv++)
        {
          //determina o centro {(h,v)} do pixel sistema de coordenadas da câmera:
          double h = (ih - ((double)nh)/2 + 0.5)*passo;
          double v = (((double)nv)/2 - iv + 0.5)*passo;

          //calcula a intersecção {p} do raio com o plano da imagem, no espaço:
          Ponto3d p = ctr;
          r3_mix_in(+h, &r, &p);
          r3_mix_in(-v, &s, &p);

          //determina a direcao do raio:
          Vetor3d dir = DirecaoDePontoParaPonto(&(cam->olho), &p);
          
          //traça o raio e armazena o resultado (seja qual for):
          Sitio sR;
          double dR;
          SitioDeInterseccaoRaioObjetos(obs, &(cam->olho), &dir, &sR, &dR);
          assert(ntot < sts.ne);
          assert(ntot%nh == ih);
          assert(ntot/nh == iv);
          sts.e[ntot] = sR;
          ntot++; 
        }
    assert(ntot == sts.ne);
    printf("gerados %d sítios de pixéis\n", ntot);

    return sts;
  }

Sitio GeraSitioSecundario(ListaDeObjetos *obs, Sitio *sO)
  {
    //perturba o sítio de origem ligeiramente na direção da normal:
    Ponto3d origem;
    r3_mix(1.0, &(sO->posicao), EPSILON_DE_ORIGEM_DO_RAIO, &(sO->normal), &origem);

    //gera uma direcao aleatoria no hemisfério positivo de {sO},
    //com distribuição enviesada para a normal de {sO}:
    Vetor3d dir;
    r3_throw_dir(&dir);
    r3_add(&(sO->normal), &dir, &dir);
    r3_mix_in(1.0e-12, &(sO->normal), &dir); //paranoia, para evitar direção nula.
    r3_dir(&dir, &dir);
    
    assert(r3_dot(&(sO->normal), &dir) >= 0.0);
    
    //lança o raio de {sO} na direção {dir}: 
    Sitio sR;
    double dR;
    SitioDeInterseccaoRaioObjetos(obs, &origem, &dir, &sR, &dR);

    if ((! isfinite(dR)) || (r3_dist(&(sR.posicao), &(sO->posicao)) <= TOLERANCIA_VISIB))
      { //sítio no infinito ou auto-intersecção, descarta:
        sR.iob = -1;
        sR.posicao = (Ponto3d){{ 0,0,0 }};
        sR.normal = (Vetor3d){{ 0,0,0 }};
      }
    
    return sR;
  }