Subversion Repositories gelsvn

Rev

Rev 346 | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 346 Rev 348
Line 1... Line 1...
1 
#ifndef __PATHTRACER_CORE__HPP__
 1 
#ifndef __PATHTRACER_CORE__HPP__
2 
#define __PATHTRACER_CORE__HPP__
 2 
#define __PATHTRACER_CORE__HPP__
3 
 
 3 
 
4 
#include "CGLA/Vec2f.h"
 4 
#include "CGLA/Vec2f.h"
5 
#include "CGLA/Vec3f.h"
 5 
#include "CGLA/Vec3f.h"
6 
 
 6 
 
7 
#include <limits>
 7 
#include <limits>
8 
 
 8 
 
9 
struct ray
 9 
struct ray
10 
{
 10 
{
11 
	ray(void) : origin(0.f), direction(0.f), distance(0.f), depth(0) {}
 11 
    ray(void) : origin(0.f), direction(0.f), distance(0.f), depth(0) {}
12 
 
 12 
 
13 
	CGLA::Vec3f origin;
 13 
    CGLA::Vec3f origin;
14 
	CGLA::Vec3f direction;
 14 
    CGLA::Vec3f direction;
15 
 
 15 
 
16 
	float distance;
 16 
    float distance;
17 
	int depth;
 17 
    int depth;
18 
};
 18 
};
19 
 
 19 
 
20 
struct hit_info
 20 
struct hit_info
21 
{
 21 
{
22 
	hit_info(void) :
22 
    hit_info(void) :
23 
		distance(0.f),
 23 
        distance(0.f),
24 
		position(0.f), geometric_normal(0.f), shading_normal(0.f),
 24 
        position(0.f), geometric_normal(0.f), shading_normal(0.f),
25 
		inside(false),
 25 
        inside(false),
26 
		emitted(0.f),
 26 
        emitted(0.f),
27 
		diffuse(0.f), glossy(0.f), shininess(0.f),
 27 
        diffuse(0.f), glossy(0.f), shininess(0.f),
28 
		reflection(0.f), refraction(0.f),
 28 
        reflection(0.f), refraction(0.f),
29 
		ior(0.f), extinction(0.f) {}
 29 
        ior(0.f), extinction(0.f) {}
30 
 
 30 
 
31 
	float distance;
 31 
    float distance;
32 
	CGLA::Vec3f position;
 32 
    CGLA::Vec3f position;
33 
	CGLA::Vec3f geometric_normal;
 33 
    CGLA::Vec3f geometric_normal;
34 
	CGLA::Vec3f shading_normal;
 34 
    CGLA::Vec3f shading_normal;
35
35 
 
36 
	//The tangent, bitangent, and shading_normal form an orthonormal basis
 36 
    //The tangent, bitangent, and shading_normal form an orthonormal basis
37 
	CGLA::Vec3f tangent;
 37 
    CGLA::Vec3f tangent;
38 
	CGLA::Vec3f bitangent;
 38 
    CGLA::Vec3f bitangent;
39 
 
 39 
 
40 
	//texture coordinates
 40 
    //texture coordinates
41 
	CGLA::Vec2f texcoords;
 41 
    CGLA::Vec2f texcoords;
42 
 
 42 
 
43 
	//flag to indicate if the inside of the object (wrt the normal) was hit
 43 
    //flag to indicate if the inside of the object (wrt the normal) was hit
44 
	bool inside;
 44 
    bool inside;
45 
 
 45 
 
46 
	//emitted radiance
 46 
    //emitted radiance
47 
	CGLA::Vec3f emitted;
 47 
    CGLA::Vec3f emitted;
48 
 
 48 
 
49 
	//diffuse reflectivity (rho)
 49 
    //diffuse reflectivity (rho)
50 
	CGLA::Vec3f diffuse;
 50 
    CGLA::Vec3f diffuse;
51
51 
 
52 
	//parameters for the modified Phong reflectance model
 52 
    //parameters for the modified Phong reflectance model
53 
	CGLA::Vec3f glossy;
 53 
    CGLA::Vec3f glossy;
54 
	float shininess;
 54 
    float shininess;
55 
 
 55 
 
56 
	//parameters for perfect specular reflection/refraction
 56 
    //parameters for perfect specular reflection/refraction
57 
	CGLA::Vec3f reflection;
 57 
    CGLA::Vec3f reflection;
58 
	CGLA::Vec3f refraction;
 58 
    CGLA::Vec3f refraction;
59 
 
 59 
 
60 
	//index of refraction + complex part
 60 
    //index of refraction + complex part
61 
	float ior;
 61 
    float ior;
62 
	float extinction;
 62 
    float extinction;
63 
};
 63 
};
64 
 
 64 
 
65 
//RGB to intensity
 65 
//RGB to intensity
66 
inline float intensity(const CGLA::Vec3f& v)
 66 
inline float intensity(const CGLA::Vec3f& v)
67 
{
 67 
{
68 
	return v[0] * 0.27f + v[1] * 0.67f + v[2] * 0.06f;
 68 
    return v[0] * 0.27f + v[1] * 0.67f + v[2] * 0.06f;
69 
}
 69 
}
70 
 
 70 
 
71 
//standard fresnel formula (see eg. Jensen p. 23)
 71 
//standard fresnel formula (see eg. Jensen p. 23)
72 
inline float fresnel_dielectric(float cosi, float cost, float eta)
 72 
inline float fresnel_dielectric(float cosi, float cost, float eta)
73 
{
 73 
{
74 
	assert(cosi>0.f && cost>0.f);
 74 
    assert(cosi>0.f && cost>0.f);
75 
 
 75 
 
76 
	float etai = 1.f;
 76 
    float etai = 1.f;
77 
	float etat = eta;
 77 
    float etat = eta;
78 
 
 78 
 
79 
	float Fs = (etai*cosi - etat*cost) / (etai*cosi + etat*cost);
 79 
    float Fs = (etai*cosi - etat*cost) / (etai*cosi + etat*cost);
80 
	float Fp = (etat*cosi - etai*cost) / (etat*cosi + etai*cost);
 80 
    float Fp = (etat*cosi - etai*cost) / (etat*cosi + etai*cost);
81
81 
 
82 
	float F = 0.5f * (Fs*Fs + Fp*Fp);
 82 
    float F = 0.5f * (Fs*Fs + Fp*Fp);
83 
	assert(F>=0.f && F<=1.f);
 83 
    assert(F>=0.f && F<=1.f);
84 
	return F;
84 
    return F;
85 
}
 85 
}
86 
 
 86 
 
87 
//standard fresnel formula (see eg. Pharr et all p. 422)
 87 
//standard fresnel formula (see eg. Pharr et all p. 422)
88 
inline float fresnel_conductor(float cosi, float eta, float extinction)
 88 
inline float fresnel_conductor(float cosi, float eta, float extinction)
89 
{
 89 
{
90 
	assert(cosi >= 0.f);
 90 
    assert(cosi >= 0.f);
91 
 
 91 
 
92 
	float z = eta*eta + extinction*extinction;
 92 
    float z = eta*eta + extinction*extinction;
93 
 
 93 
 
94 
	float Fp = (z * cosi*cosi - 2.f * eta * cosi + 1.f)	/
 94 
    float Fp = (z * cosi*cosi - 2.f * eta * cosi + 1.f) /
95 
		(z + 2.f * eta * cosi + 1.f);
 95 
        (z + 2.f * eta * cosi + 1.f);
96 
 
 96 
 
97 
	float Fs = (z - 2.f * eta * cosi + cosi * cosi) /
 97 
    float Fs = (z - 2.f * eta * cosi + cosi * cosi) /
98 
		(z + 2.f * eta * cosi + cosi * cosi);
 98 
        (z + 2.f * eta * cosi + cosi * cosi);
99 
 
 99 
 
100 
	float F = (Fp + Fs) * 0.5f;
 100 
    float F = (Fp + Fs) * 0.5f;
101 
	assert(F>=0.f && F<=1.f);
 101 
    assert(F>=0.f && F<=1.f);
102 
	return F;
 102 
    return F;
103 
}
 103 
}
104 
 
 104 
 
105 
//reflect vector around normal. n should be normalized.
 105 
//reflect vector around normal. n should be normalized.
106 
inline CGLA::Vec3f reflect(const CGLA::Vec3f& n, const CGLA::Vec3f& r)
 106 
inline CGLA::Vec3f reflect(const CGLA::Vec3f& n, const CGLA::Vec3f& r)
107 
{
 107 
{
108 
	return 2.f * n * dot(n, r) - r;
 108 
    return 2.f * n * dot(n, r) - r;
109 
}
 109 
}
110 
 
 110 
 
111 
//eta = eta_from / eta_to
 111 
//eta = eta_from / eta_to
112 
inline bool refract(const CGLA::Vec3f & n, const CGLA::Vec3f & i,
112 
inline bool refract(const CGLA::Vec3f & n, const CGLA::Vec3f & i,
113 
					float eta, CGLA::Vec3f & t)
 113 
                    float eta, CGLA::Vec3f & t)
114 
{
 114 
{
115 
	float c1 = dot(i,n);
 115 
    float c1 = dot(i,n);
116 
	float c2 = 1.f - eta*eta * (1.f - c1*c1);
 116 
    float c2 = 1.f - eta*eta * (1.f - c1*c1);
117 
 
 117 
 
118 
	if (c2 < 0.f)
 118 
    if (c2 < 0.f)
119 
		return false;
 119 
        return false;
120 
 
 120 
 
121 
	float c3 = std::sqrt(c2);
 121 
    float c3 = std::sqrt(c2);
122 
	t = -eta*i + (eta*c1 - c3) * n;
 122 
    t = -eta*i + (eta*c1 - c3) * n;
123 
 
 123 
 
124 
	return true;
 124 
    return true;
125 
}
 125 
}
126 
 
 126 
 
127 
//helper function, clamps v to range [i; a]
 127 
//helper function, clamps v to range [i; a]
128 
template <class T>
 128 
template <class T>
129 
T clamp(const T& v, const T& i, const T& a)
 129 
T clamp(const T& v, const T& i, const T& a)
130 
{
 130 
{
131 
	return std::min(std::max(v, i), a);
 131 
    return std::min(std::max(v, i), a);
132 
}
 132 
}
133 
 
 133 
 
134 
//mersenne twister
 134 
//mersenne twister
135 
extern "C" double genrand_real2(void);
 135 
double genrand_real2(void);
136 
inline float random(void)
 136 
inline float mt_random(void)
137 
{
 137 
{
138 
	return float(genrand_real2());
 138 
    return float(genrand_real2());
139 
}
 139 
}
140 
 
 140 
 
141 
//evaluation of the diffuse brdf
 141 
//evaluation of the diffuse brdf
142 
inline CGLA::Vec3f lambertian_brdf(const hit_info& hi,
 142 
inline CGLA::Vec3f lambertian_brdf(const hit_info& hi,
143 
								   const CGLA::Vec3f& wi,
 143 
                                   const CGLA::Vec3f& wi,
144 
								   const CGLA::Vec3f& wo)
 144 
                                   const CGLA::Vec3f& wo)
145 
{
 145 
{
146 
	return hi.diffuse / float(M_PI);
 146 
    return hi.diffuse / float(M_PI);
147 
}
 147 
}
148 
 
 148 
 
149 
//sampling of diffuse brdf, wi = sampled direction, should return probability
 149 
//sampling of diffuse brdf, wi = sampled direction, should return probability
150 
inline float sample_lambertian(const hit_info& hi,
 150 
inline float sample_lambertian(const hit_info& hi,
151 
							const CGLA::Vec3f& wo, CGLA::Vec3f& wi)
 151 
                            const CGLA::Vec3f& wo, CGLA::Vec3f& wi)
152 
{
 152 
{
153 
	CGLA::Vec3f x, y, z = hi.shading_normal;
 153 
    CGLA::Vec3f x, y, z = hi.shading_normal;
154 
	orthogonal(z, x, y);
 154 
    orthogonal(z, x, y);
155 
	float e0 = random();
 155 
    float e0 = mt_random();
156 
	float e1 = random();
 156 
    float e1 = mt_random();
157 
	float cost = std::sqrt(e0);
 157 
    float cost = std::sqrt(e0);
158 
	float sint = std::sqrt(1.f - e0);
 158 
    float sint = std::sqrt(1.f - e0);
159 
	float cosp = std::cos(e1 * 2.f * float(M_PI));
 159 
    float cosp = std::cos(e1 * 2.f * float(M_PI));
160 
	float sinp = std::sin(e1 * 2.f * float(M_PI));
 160 
    float sinp = std::sin(e1 * 2.f * float(M_PI));
161 
	CGLA::Vec3f dir(sint*cosp, sint*sinp, cost);
 161 
    CGLA::Vec3f dir(sint*cosp, sint*sinp, cost);
162 
 
 162 
 
163 
	wi = dir(0) * x + dir(1) * y + dir(2) * z;
 163 
    wi = dir(0) * x + dir(1) * y + dir(2) * z;
164 
	return 1.f / float(M_PI) * cost;
 164 
    return 1.f / float(M_PI) * cost;
165 
}
 165 
}
166 
 
 166 
 
167 
//evaluation of the phong brdf
167 
//evaluation of the phong brdf
168 
inline CGLA::Vec3f phong_brdf(const hit_info& hi,
 168 
inline CGLA::Vec3f phong_brdf(const hit_info& hi,
169 
							   const CGLA::Vec3f& wi,const CGLA::Vec3f& wo)
 169 
                               const CGLA::Vec3f& wi,const CGLA::Vec3f& wo)
170 
{
 170 
{
171 
	float dwi = dot(wi, hi.geometric_normal);
 171 
    float dwi = dot(wi, hi.geometric_normal);
172 
	float dwo = dot(wo, hi.geometric_normal);
 172 
    float dwo = dot(wo, hi.geometric_normal);
173 
	float same_hemisphere = dwi * dwo > 0.f;
 173 
    float same_hemisphere = dwi * dwo > 0.f;
174 
 
 174 
 
175 
	if (same_hemisphere)
 175 
    if (same_hemisphere)
176 
	{
 176 
    {
177 
		//glossy contribution
 177 
        //glossy contribution
178 
		if (intensity(hi.glossy) > 0.f)
 178 
        if (intensity(hi.glossy) > 0.f)
179 
		{
 179 
        {
180 
			float cost = std::max(dot(reflect(hi.shading_normal, wi), wo),
 180 
            float cost = std::max(dot(reflect(hi.shading_normal, wi), wo),
181 
				0.f);
 181 
                0.f);
182 
			float cosn = std::pow(cost, hi.shininess);
 182 
            float cosn = std::pow(cost, hi.shininess);
183 
			return hi.glossy * (hi.shininess+2.f)/(2.f * float(M_PI))*cosn;
 183 
            return hi.glossy * (hi.shininess+2.f)/(2.f * float(M_PI))*cosn;
184 
		}
 184 
        }
185 
	}
 185 
    }
186 
	return CGLA::Vec3f(0.f);
 186 
    return CGLA::Vec3f(0.f);
187 
}
 187 
}
188 
 
 188 
 
189 
//sampling of phong brdf, wi = sampled direction, should return probability
 189 
//sampling of phong brdf, wi = sampled direction, should return probability
190 
inline float sample_phong(const hit_info& hi,
 190 
inline float sample_phong(const hit_info& hi,
191 
						  const CGLA::Vec3f& wo, CGLA::Vec3f& wi)
 191 
                          const CGLA::Vec3f& wo, CGLA::Vec3f& wi)
192 
{
 192 
{
193 
	float e0 = random();
 193 
    float e0 = mt_random();
194 
	float e1 = random();
 194 
    float e1 = mt_random();
195
195 
 
196 
	float cost = std::pow(e0, 1.f/(hi.shininess + 1.f));
 196 
    float cost = std::pow(e0, 1.f/(hi.shininess + 1.f));
197 
	float sint = std::sqrt(1.f - std::pow(e0, 2.f/(hi.shininess + 1.f)));
 197 
    float sint = std::sqrt(1.f - std::pow(e0, 2.f/(hi.shininess + 1.f)));
198 
	float cosp = std::cos(e1 * 2.f * float(M_PI));
 198 
    float cosp = std::cos(e1 * 2.f * float(M_PI));
199 
	float sinp = std::sin(e1 * 2.f * float(M_PI));
 199 
    float sinp = std::sin(e1 * 2.f * float(M_PI));
200
200 
 
201 
	CGLA::Vec3f dir(sint*cosp, sint*sinp, cost);
 201 
    CGLA::Vec3f dir(sint*cosp, sint*sinp, cost);
202 
 
 202 
 
203 
	CGLA::Vec3f x, y, z = reflect(hi.shading_normal, wo);
 203 
    CGLA::Vec3f x, y, z = reflect(hi.shading_normal, wo);
204 
	orthogonal(z, x, y);
 204 
    orthogonal(z, x, y);
205 
 
 205 
 
206 
	float cosn = std::pow(dir(2), hi.shininess);
 206 
    float cosn = std::pow(dir(2), hi.shininess);
207 
	wi = dir(0) * x + dir(1) * y + dir(2) * z;
 207 
    wi = dir(0) * x + dir(1) * y + dir(2) * z;
208 
 
 208 
 
209 
	return (hi.shininess+1.f) / (2.f * float(M_PI)) * cosn;
 209 
    return (hi.shininess+1.f) / (2.f * float(M_PI)) * cosn;
210 
}
 210 
}
211 
 
 211 
 
212 
//convenience function to evaluate the non-specular parts of the bsdf
 212 
//convenience function to evaluate the non-specular parts of the bsdf
213 
inline CGLA::Vec3f bsdf_evaluate(const hit_info& hi,
213 
inline CGLA::Vec3f bsdf_evaluate(const hit_info& hi,
214 
								 const CGLA::Vec3f& wi, const CGLA::Vec3f& wo)
 214 
                                 const CGLA::Vec3f& wi, const CGLA::Vec3f& wo)
215 
{
 215 
{
216 
	CGLA::Vec3f fs(0.f);
 216 
    CGLA::Vec3f fs(0.f);
217 
 
 217 
 
218 
	fs += lambertian_brdf(hi, wi, wo);
 218 
    fs += lambertian_brdf(hi, wi, wo);
219 
	fs += phong_brdf(hi, wi, wo);
 219 
    fs += phong_brdf(hi, wi, wo);
220 
 
 220 
 
221 
	return fs;
 221 
    return fs;
222 
}
 222 
}
223 
 
 223 
 
224 
//epsilon for shadow testing
 224 
//epsilon for shadow testing
225 
static const float epsilon = 1e-5f;
 225 
static const float epsilon = 1e-5f;
226 
 
 226