#include "tile_renderer.h"

#include "renderer.h"
#include "../noise/Simplex.h"

DT_FP mod(DT_FP a, DT_FP b) {
	DT_FP result = fmod(a,b);
	if (result<(DT_FP)(0.0)) result += b;
	return result;
}

void add_hsv_to_rgb(DT_COL hsv[3], DT_COL rgb[3]) {
	DT_COL chroma = hsv[2] * hsv[1];
	DT_COL hprime = hsv[0] / 60.0;
	DT_COL x = chroma  *  ( (DT_COL)(1.0) - fabs(fmod(hprime,(DT_COL)(2.0))-(DT_COL)(1.0)) );
	if      ( hprime>=(DT_COL)(0.0) && hprime<(DT_COL)(1.0) ) { rgb[0] +=        chroma; rgb[1] +=             x; rgb[2] += (DT_COL)(0.0); }
	else if ( hprime>=(DT_COL)(1.0) && hprime<(DT_COL)(2.0) ) { rgb[0] +=             x; rgb[1] +=        chroma; rgb[2] += (DT_COL)(0.0); }
	else if ( hprime>=(DT_COL)(2.0) && hprime<(DT_COL)(3.0) ) { rgb[0] += (DT_COL)(0.0); rgb[1] +=        chroma; rgb[2] +=             x; }
	else if ( hprime>=(DT_COL)(3.0) && hprime<(DT_COL)(4.0) ) { rgb[0] += (DT_COL)(0.0); rgb[1] +=             x; rgb[2] +=        chroma; }
	else if ( hprime>=(DT_COL)(4.0) && hprime<(DT_COL)(5.0) ) { rgb[0] +=             x; rgb[1] += (DT_COL)(0.0); rgb[2] +=        chroma; }
	else if ( hprime>=(DT_COL)(5.0) && hprime<(DT_COL)(6.0) ) { rgb[0] +=        chroma; rgb[1] += (DT_COL)(0.0); rgb[2] +=             x; }
	else {
		throw glLib::error("Something went wrong in converting from HSV!");
	}
}
void add_color_heatmap(DT_COL color[3], DT_COL intensity) {
	intensity = glLib::Math::Misc::clamp(intensity,(DT_COL)(0.0),(DT_COL)(1.0));
	DT_COL hsv[3] = {((DT_COL)(1.0)-intensity)*(DT_COL)(240.0),(DT_COL)(1.0),(DT_COL)(1.0)};
	add_hsv_to_rgb(hsv,color);
}
void map_vec_to_rgb(DT_FP vec[3], DT_COL color[3]) {
	#ifdef OUTPUT_ABS
	color[0] = fabs(vec[0]);
	color[1] = fabs(vec[1]);
	color[2] = fabs(vec[2]);
	#else
	color[0] = (vec[0]+1.0)*0.5;
	color[1] = (vec[1]+1.0)*0.5;
	color[2] = (vec[2]+1.0)*0.5;
	#endif
}

DT_FP fresnel_rs(DT_FP n1, DT_FP n2, DT_FP theta_i, DT_FP theta_t) {
	DT_FP a = n1*cos(theta_i), b = n2*cos(theta_t);
	DT_FP term = (a-b) / (a+b);
	return term * term;
}
DT_FP fresnel_rp(DT_FP n1, DT_FP n2, DT_FP theta_i, DT_FP theta_t) {
	DT_FP a = n1*cos(theta_t), b = n2*cos(theta_i);
	DT_FP term = (a-b) / (a+b);
	return term * term;
}
inline DT_FP averaged_fresnel_r(DT_FP n1, DT_FP n2, DT_FP theta_i, DT_FP theta_t) {
	return ( fresnel_rs(n1,n2,theta_i,theta_t) + fresnel_rp(n1,n2,theta_i,theta_t) )  /  (DT_FP)(2.0);
}
void averaged_fresnel(DT_FP n1,DT_FP n2,DT_FP theta_i,DT_FP theta_t, DT_FP coefficient_reflection[],DT_FP coefficient_transmission[]) {
	coefficient_reflection[0] = averaged_fresnel_r(n1,n2,theta_i,theta_t);
	coefficient_transmission[0] = (DT_FP)(1.0) - coefficient_reflection[0];
	for (char i=1;i<3;++i) {
		coefficient_reflection[i] = coefficient_reflection[0];
		coefficient_transmission[i] = coefficient_transmission[0];
	}
}
void averaged_fresnel_transmit(DT_FP n1,DT_FP n2,DT_FP theta_i,DT_FP theta_t, DT_FP coefficient_transmission[]) {
	coefficient_transmission[0] = (DT_FP)(1.0) - averaged_fresnel_r(n1,n2,theta_i,theta_t);
	for (char i=1;i<3;++i) {
		coefficient_transmission[i] = coefficient_transmission[0];
	}
}
void averaged_fresnel_reflect(DT_FP n1,DT_FP n2,DT_FP theta_i,DT_FP theta_t, DT_FP coefficient_reflection[]) {
	coefficient_reflection[0] = averaged_fresnel_r(n1,n2,theta_i,theta_t);
	for (char i=1;i<3;++i) {
		coefficient_reflection[i] = coefficient_reflection[0];
	}
}

void print_color(DT_COL color[3]) {
	printf("(");
	for (char i=0;i<3;++i) {
		printf("%3d",glLib::Math::Misc::rnduint(glLib::Math::Misc::clamp<DT_COL>(255.0*color[i],0,255)));
		if (i!=2) printf(", ");
	}
	printf(")");
}

DT_FP interfere_photon_amplitude(DT_FP amp1, DT_FP amp2, DT_FP lambda, DT_FP shift) {
	//Returns the new amplitude.  Note that the shift might change also.
	return fabs(  amp1 + amp2 * cos( glLib::Math::Misc::two_pi * shift/lambda )  );
}

void thickness_offset(DT_FP& film_thickness, DT_FP position_intersection[3]) {
	const static int octaves = 5;
	const static double persistence = 0.5;
	const static double scale = 0.5;

	double sample000 = noise_simplex(
		octaves,persistence,scale, position_intersection[0],position_intersection[1],position_intersection[2]
	);
#if 0
	double sample100 = noise_simplex(
		octaves,persistence,scale, position_intersection[0]+EPSILON,position_intersection[1],position_intersection[2]
	);
	double sample010 = noise_simplex(
		octaves,persistence,scale, position_intersection[0],position_intersection[1]+EPSILON,position_intersection[2]
	);
	double sample001 = noise_simplex(
		octaves,persistence,scale, position_intersection[0],position_intersection[1],position_intersection[2]+EPSILON
	);

	//Calculate curl
	double dx=sample100-sample000, dy=sample010-sample000, dz=sample001-sample000;
	double curl[3] = {dy-dz,dz-dx,dx-dy};
	glLib::Math::Vec::normalize3(curl);

	//Offset at most 100% of the film's original thickness
	film_thickness += 1.00*film_thickness*curl[2];
#else
	//Offset at most 100% of the film's original thickness
	film_thickness += 1.0*film_thickness*sample000;
#endif
}
//Fully refractive
DT_FP thin_film_reflect_refract(DT_FP n1, DT_FP n2, DT_FP theta_1, DT_FP film_index, DT_FP film_thickness, DT_FP position_intersection[3], DT_FP coefficient_reflection[3], DT_FP coefficient_transmission[3]) {
	//Reflection 1st interaction: coeff1
	//Refraction 1st interaction: coeff2
	//Reflection 2nd interaction: coeff3
	//Refraction 2nd interaction: coeff4
	//Reflection 3rd interaction: coeff5 (ignored, so doesn't exist)
	//Refraction 3rd interaction: coeff6
	DT_FP coeff1[3], coeff2[3], coeff3[3], coeff4[3], coeff6[3];

	//1st interaction (incident hitting film)
	DT_FP theta_2 = asin( (n1/film_index) * sin(theta_1) ); //theta2 = arcsin((n1/n2)*sin(theta1))
	averaged_fresnel( n1,film_index, theta_1,theta_2, coeff1,coeff2 );

	//2nd interaction (transmitted hitting opposite side of film)
	DT_FP theta_3 = asin( (film_index/n2) * sin(theta_2) ); //theta2 = arcsin((n1/n2)*sin(theta1))
	averaged_fresnel( film_index,n2, theta_2,theta_3, coeff3,coeff4 );
	//	Must be scaled by transmitted amount.
	glLib::Math::Vec::multiply3(coeff2, coeff3,coeff3);
	glLib::Math::Vec::multiply3(coeff2, coeff4,coeff4);

	//3rd interaction (reflected off opposite side, refracting back out to n1)
	averaged_fresnel_transmit( film_index,n1, theta_2,theta_1, coeff6 );
	//	Must be scaled by amount reflected off opposite side.
	glLib::Math::Vec::multiply3(coeff3, coeff6,coeff6);

	thickness_offset(film_thickness,position_intersection);
	DT_FP path_length = 2.0*film_thickness / cos(theta_2);

	//Reflection coefficient has thin film interference (basically, the whole point of this)
	coefficient_reflection[0] = interfere_photon_amplitude( coeff1[0],coeff6[0], 680, path_length );
	coefficient_reflection[1] = interfere_photon_amplitude( coeff1[1],coeff6[1], 545, path_length );
	coefficient_reflection[2] = interfere_photon_amplitude( coeff1[2],coeff6[2], 465, path_length );
	for (char i=0;i<3;++i) {
		coefficient_transmission[i] = (DT_FP)(1.0) - coefficient_reflection[i];
	}

	return path_length;
}
//Assumes completely reflective underlying substrate
DT_FP thin_film_reflect        (DT_FP n1,           DT_FP theta_1, DT_FP film_index, DT_FP film_thickness, DT_FP position_intersection[3], DT_FP coefficient_reflection[3]                                   ) {
	//Reflection 1st interaction: coeff1
	//Refraction 1st interaction: coeff2  \
	//Reflection 2nd interaction: coeff3  _\same, so named coeff23
	//Refraction 2nd interaction: coeff4 (ignored, so doesn't exist)
	//Reflection 3rd interaction: coeff5 (ignored, so doesn't exist)
	//Refraction 3rd interaction: coeff6
	DT_FP coeff1[3], coeff23[3], coeff6[3];

	//1st interaction (incident hitting film)
	DT_FP theta_2 = asin( (n1/film_index) * sin(theta_1) ); //theta2 = arcsin((n1/n2)*sin(theta1))
	averaged_fresnel( n1,film_index, theta_1,theta_2, coeff1,coeff23 );

	//2nd interaction (transmitted hitting opposite side of film)
	//Does nothing; total internal reflection back into film

	//3rd interaction (reflected off opposite side, refracting back out to n1)
	averaged_fresnel_transmit( film_index,n1, theta_2,theta_1, coeff6 );
	//	Must be scaled by amount reflected off opposite side.
	glLib::Math::Vec::multiply3(coeff23, coeff6,coeff6);

	thickness_offset(film_thickness,position_intersection);
	DT_FP path_length = 2.0*film_thickness / cos(theta_2);

	//Reflection coefficient has thin film interference (basically, the whole point of this)
	coefficient_reflection[0] = interfere_photon_amplitude( coeff1[0],coeff6[0], 680, path_length );
	coefficient_reflection[1] = interfere_photon_amplitude( coeff1[1],coeff6[1], 545, path_length );
	coefficient_reflection[2] = interfere_photon_amplitude( coeff1[2],coeff6[2], 465, path_length );
	
	return path_length;
}

//Numerical error caused an entry/exit to be less than an epsilon
//or normals don't match surface perfectly.  Maybe we started inside
//a surface already.  In any case, something horrible happened.
#define RETURN_ERROR { color[0]=1000.0;color[1]=0.0;color[2]=0.0;trace_info.error=true;return trace_info; }

TileRenderer* instance__renderer;
inline int wrapper__render_tile(void* arg) {
	return instance__renderer->render_tile(arg);
}

TileRenderer::TileRenderer(Renderer* renderer, int id) {
	_renderer = renderer;

	material_air = new Material("air");

	init_matrices();

	_thread = NULL;
	_thread_id = id;

	instance__renderer = this;
}
TileRenderer::~TileRenderer() {
	deinit_matrices();

	delete material_air;
}

void TileRenderer::thread_start(void) {
	_thread = SDL_CreateThread(wrapper__render_tile,&_thread_id);
}
void TileRenderer::thread_wait(void) {
	if (_thread!=NULL) SDL_WaitThread(_thread,NULL);
}

unsigned int TileRenderer::trace_get_intersect(Ray ray, RayInfo* info) {
	#ifdef USING_BVH
		return _renderer->_bvh->intersect(ray,info);
	#else
		info->reset();
		for (list<Object*>::const_iterator iter=_renderer->_scene->objects->begin(); iter!=_renderer->_scene->objects->end(); ++iter) {
			RayInfo temp = RayInfo();
			bool hit = intersect_object(iter,ray,&temp);
			if (hit) {
				if (temp.hit_distance<=0.0) continue;
				if (temp.hit_distance<info->hit_distance) {
					info->hit_object = *iter;
					info->set_hit(*iter,temp.hit_distance);
				}
			}
		}
		return _renderer->_scene->objects->size();
	#endif
}
#ifdef WHITTED
RayInfo TileRenderer::trace_shadow_ray(Ray ray) {
	RayInfo trace_info;
	trace_get_intersect(ray,&trace_info);
	return trace_info;
}
RayInfo TileRenderer::trace_ray(Ray ray_incident, DT_COL color[3], int depth, deque<Material*>* material_stack) {
	//Determine intersect
	RayInfo trace_info;
	unsigned int num_tests = trace_get_intersect(ray_incident,&trace_info);

	//Set the materials
	if (material_stack->empty()) RETURN_ERROR
	Material* material_current = material_stack->front();

	#ifdef DEBUG_CLICK
	if (_renderer->render_debug) {
		//DEBUG_CLICK_INDENT printf("-Traveling through %s\n",material_current->name->c_str());
		//DEBUG_CLICK_INDENT printf("-With absorption coefficient "); glLib::Math::Vec::print3(material_current->absorption);
		//DEBUG_CLICK_INDENT
		//if (trace_info.hit_object!=NULL) printf("-Hit object \"%s\"\n",trace_info.hit_object->name->c_str());
		//else printf("-Hit object (None)\n");
	}
	#endif

	#ifdef OUTPUT_BVH_HEATMAP
	add_color_heatmap(color,num_tests/15.0);
	return trace_info;
	#endif

	//Shade
	if (trace_info.hit) {
		#ifdef OUTPUT_T
		add_color_heatmap( (closest-OUTPUT_T_NEAR)/(OUTPUT_T_FAR-OUTPUT_T_NEAR) );
		return trace_info;
		#endif

		#ifdef ALLOW_AMBIENT
		glLib::Math::Vec::add3(trace_info.hit_material->ambient,color,color);
		#endif

		DT_FP position_intersection[3];
		glLib::Math::Vec::scale3(trace_info.hit_distance,ray_incident.get_direction(),position_intersection);
		glLib::Math::Vec::add3(position_intersection,ray_incident.origin,position_intersection);

		DT_FP vec_light[3];
		DT_FP vec_normal[3];
		DT_FP vec_reflection[3];
		DT_FP vec_refraction[3];
		DT_FP vec_source[3];

		glLib::Math::Vec::copy3(ray_incident.get_direction_negative(),vec_source);

		//Get normal
		trace_info.hit_object->set_normal(position_intersection,vec_normal);

		#ifdef OUTPUT_NORMAL
		map_vec_to_rgb(vec_normal,color);
		return trace_info;
		#endif

		#if defined ALLOW_DIFFUSE || defined ALLOW_SPECULAR// || defined ALLOW_SHADOWS
		map<string,Light*>::const_iterator iter;
		for (iter=_renderer->_lights->begin();iter!=_renderer->_lights->end();++iter) {
			Light* light = iter->second;
			glLib::Math::Vec::subtract3(light->primitive->position,position_intersection,vec_light);
			glLib::Math::Vec::normalize3(vec_light);

			DT_FP n_dot_l = glLib::Math::Vec::dot3(vec_normal,vec_light);
			if (n_dot_l<(DT_FP)(0.0)) {
				//if (double_sided) {
				//	glLibMathVec::negate3(n);
				//	n_dot_l = -n_dot_l;
				//} else {
				continue;
				//}
			}

			DT_COL light_dependent_color[3], light_dependent_scalar = (DT_FP)(1.0);

			#ifdef ALLOW_DIFFUSE
			//Diffuse
			DT_COL diffuse_color[3]; glLib::Math::Vec::copy3(trace_info.hit_material->diffuse,diffuse_color);
			SDL_Surface* texture = trace_info.hit_material->texture;
			switch (trace_info.hit_object->type) {
				case SCENE_OBJECT_TYPE_PLANE: {
					DT_FP sc = 4.0;
					if (((Plane*)(trace_info.hit_object))->normal[1]==1) {
						if ((mod(sc*position_intersection[0],2.0)>1.0)^(mod(sc*position_intersection[2],2.0)>1.0)) {
							n_dot_l *= 0.5;//*max(mod(sc*position_intersection[0],1.0),mod(sc*position_intersection[2],1.0));
						}
					}
					else {
						if ((mod(sc*position_intersection[0],2.0)>1.0)^(mod(sc*position_intersection[1],2.0)>1.0)) {
							n_dot_l *= 0.5;//9*max(mod(sc*position_intersection[0],1.0),mod(sc*position_intersection[1],1.0));
						}
					}
					break;
				}
				case SCENE_OBJECT_TYPE_SPHERE: {
					if (texture!=NULL) {
						DT_FP uv[2] = {0.5*texture->w,0.5*texture->h};

						DT_FP angle_x = atan2(vec_normal[2],vec_normal[0]);
						DT_FP angle_y = atan2(vec_normal[1],sqrt(vec_normal[0]*vec_normal[0]+vec_normal[2]*vec_normal[2]));
						uv[0] = texture->w * (angle_x + glLib::Math::Misc::pi    )/(2.0*glLib::Math::Misc::pi);
						uv[1] = texture->h * (angle_y + glLib::Math::Misc::pi*0.5)/     glLib::Math::Misc::pi ;
						DT_COL texture_color[3];
						sample_default(uv[0],uv[1],texture,texture_color);
						glLib::Math::Vec::multiply3(texture_color,diffuse_color,diffuse_color);
					}
					break;
				}
				case SCENE_OBJECT_TYPE_TRIANGLE: {
					if (texture!=NULL) {
						DT_FP uv[2] = {0.5*texture->w,0.5*texture->h};

						Triangle* t = (Triangle*)(trace_info.hit_object);
						t->set_uv(position_intersection,uv);
						DT_COL texture_color[3];
						sample_default(uv[0]*texture->w,uv[1]*texture->h,texture,texture_color);
						glLib::Math::Vec::multiply3(texture_color,diffuse_color,diffuse_color);
					}
					break;
				}
				default: throw glLib::error("Unrecognized type!");
			}

			#ifdef OUTPUT_UV_DIFFUSE
			color[0] = uv[0]/texture->w;
			color[1] = uv[1]/texture->h;
			color[2] = 0.0;
			return trace_info;
			#endif

			DT_COL scaled_diffuse[3]; glLib::Math::Vec::scale3(n_dot_l,diffuse_color,scaled_diffuse);
			glLib::Math::Vec::copy3(scaled_diffuse,light_dependent_color);
			#endif

			#ifdef ALLOW_SPECULAR
			//Specular
			//	Not using glLibMathVec::reflect, because can take advantage of the fact that dot(n,l) has already been computed.
			DT_FP temp[3];
			glLib::Math::Vec::scale3((DT_FP)(2.0)*n_dot_l,vec_normal,temp);
			glLib::Math::Vec::subtract3(temp,vec_light,vec_reflection);

			DT_FP spec = glLib::Math::Vec::dot3(vec_reflection,vec_source);
			spec = (spec<(DT_FP)(0.0)) ? (DT_FP)(0.0) : spec; //when raising small negative numbers to large powers, sometimes it becomes positive.
			//spec = pow(spec,material_hit->shininess);

			//21.85007974481659% faster
			int pow_index = glLib::Math::Misc::clamp( (int)(spec*(NUM_POW_TABLE-1)), 0, NUM_POW_TABLE-2); //inaccuracies?
			spec = (trace_info.hit_material->pow_table[pow_index]+trace_info.hit_material->pow_table[pow_index+1])/2.0;

			DT_COL scaled_specular[3]; glLib::Math::Vec::scale3(spec,trace_info.hit_material->specular,scaled_specular);
			glLib::Math::Vec::add3(scaled_specular,light_dependent_color,light_dependent_color);
			#endif

			#ifdef ALLOW_SHADOWS
			//Shadow
			map<Light*,int>::const_iterator location = trace_info.hit_material->to_shadow->find(light);
			if (location!=trace_info.hit_material->to_shadow->end()) {
				light_dependent_scalar = 0.0;
				for (int i=0;i<location->second;++i) {
					DT_FP light_point[3];
					light->get_random_point(light_point);

					DT_FP shadow_ray_direction[3];
					glLib::Math::Vec::subtract3(light_point,position_intersection,shadow_ray_direction);

					DT_FP distance_to_light_point = glLib::Math::Vec::length3(shadow_ray_direction);
					glLib::Math::Vec::scale3( (DT_FP)(1.0)/distance_to_light_point, shadow_ray_direction );

					DT_FP offset_position[3];
					glLib::Math::Vec::scaleadd3(EPSILON,vec_normal, position_intersection, offset_position);
					Ray shadow_ray(offset_position,shadow_ray_direction);
					RayInfo shadow_info = trace_shadow_ray(shadow_ray);
					if (  (!shadow_info.hit)  ||  ((shadow_info.hit)&&(shadow_info.hit_distance>distance_to_light_point)) ) {
						light_dependent_scalar += 1.0;
					}
				}
				light_dependent_scalar /= (DT_COL)(location->second);
			}
			#endif

			//Scale light dependent color
			glLib::Math::Vec::scale3(light_dependent_scalar,light_dependent_color);

			//Add diffuse and specular
			glLib::Math::Vec::multiply3(light_dependent_color,light->color,light_dependent_color);
			glLib::Math::Vec::add3(light_dependent_color,color,color);
		}
		#endif

		//Handle refractions and reflections.
		if (depth<MAX_DEPTH) {
			//Whether we're entering or leaving a material
			char entering_leaving = 1;
			//If the normal is facing away, then we're inside an object.  Flip the normal around,
			//so that we can use it, but make a note of it with entering_leaving.
			if (glLib::Math::Vec::dot3(vec_source,vec_normal)<0.0) {
				//Change the normal so as to be facing toward the incident vector.
				glLib::Math::Vec::negate3(vec_normal);
				entering_leaving = 2;
			}

			//Fresnel coefficients: n1 = 1.0, n2 = infinity <=> opaque object
			DT_FP   coefficient_reflection[3];
			DT_FP coefficient_transmission[3];

			deque<Material*> material_stack_refract_ray = *material_stack; //copies deque values
			deque<Material*> material_stack_reflect_ray = *material_stack; //copies deque values

			Material* material_new = NULL;
			if (entering_leaving==1) { //entering something
				material_new = trace_info.hit_material;
				material_stack_refract_ray.push_front(material_new);
			}
			else { //leaving something
				material_stack_refract_ray.pop_front();
				if (material_stack_refract_ray.empty()) RETURN_ERROR
				material_new = material_stack_refract_ray.front();
			}
			DT_FP  n1 = (DT_FP)(material_current->refractive_index);
			DT_FP  n2 = (DT_FP)(    material_new->refractive_index);
			DT_FP eta = n1 / n2;

			//Variables requisite for calculating Fresnel coefficients.
			DT_FP cos_theta_i = glLib::Math::Vec::dot3(vec_source,vec_normal);
			//	1, 2, 3 => incident, refracted1, refracted2
			DT_FP theta_1 = acos(cos_theta_i);//n1*sin(theta1) = n2*sin(theta2)

			#define film_index trace_info.hit_material->film_refractive_index
			#define film_thickness trace_info.hit_material->film_thickness

			bool must_calculate_fresnel = true;

			#ifdef ALLOW_REFRACTIONS
			if (trace_info.hit_material->refract) {
				//Total internal reflection; leave it for standard reflection.
				if (!glLib::Math::Vec::refract3( vec_source, vec_normal, vec_refraction, eta )) {
					#ifdef DEBUG_CLICK
					if (_renderer->render_debug){
						DEBUG_CLICK_INDENT printf("-Total internal reflection!\n");
						DEBUG_CLICK_INDENT printf("Not tracing refraction ray!\n");
					}
					#endif
					//color[0] = 0.0; color[1] = 0.0; color[2] = 10000.0;
					//return trace_info;
					// printf("Total internal reflection!  Incident %f\n",glLib::Math::Misc::degrees(acos(glLib::Math::Vec::dot3(vec_source,vec_normal))));
				}
				//Refraction was successful!  Reflection handles other ray.
				else {
					//Set the Fresnel Coefficients, attenuating by thin films, if necessary (and allowed).
					#ifdef ALLOW_FILMS
					if (trace_info.hit_material->film) {
						thin_film_reflect_refract( n1,n2, theta_1, film_index,film_thickness, position_intersection, coefficient_reflection,coefficient_transmission );
						//TODO:
						#ifdef OUTPUT_FILM_PATH_LENGTH
						add_color_heatmap(color, path_length/OUTPUT_FILM_PATH_LENGTH_SCALE); return trace_info;
						#endif
					} else {
					#endif
						//Should never give a domain error; should fail above as total internal reflection
						DT_FP theta_3 = asin( eta*sin(theta_1) ); //theta2 = arcsin((n1/n2)*sin(theta1))
						averaged_fresnel( n1,n2, theta_1,theta_3, coefficient_reflection,coefficient_transmission );
					#ifdef ALLOW_FILMS
					}
					#endif
					must_calculate_fresnel = false;

					#ifdef DEBUG_CLICK
					if (_renderer->render_debug) {
						//DEBUG_CLICK_INDENT printf("-Indices of refraction: n1=%1.4f; n2=%1.4f\n",n1,n2);
					}
					#endif

					//color[0] = 0.0; color[1] = 10000.0; color[2] = 0.0;
					//return trace_info;

					#ifdef OUTPUT_REFRACTION
					//if (depth==OUTPUT_REFRACTION_DEPTH) {
					if (depth==0) {
						map_vec_to_rgb(vec_refraction,color);
						return trace_info;
					}
					#endif

					//Sample a refraction direction
					DT_COL color_sample_refraction[3] = {0.0,0.0,0.0};

					DT_FP offset_position[3];
					glLib::Math::Vec::scaleadd3(-EPSILON,vec_normal, position_intersection, offset_position);

					/*if (trace_info.hit_material->samples_refraction>1) {
						DT_FP world_y[3] = {0.0,1.0,0.0};
						DT_FP tangent[3]; glLib::Math::Vec::cross(world_y,vec_refraction,tangent); glLib::Math::Vec::normalize3(tangent);
						DT_FP distribution_angular_offset = acos(glLib::Math::Vec::dot3(world_y,vec_refraction));

						glLib::Math::Matrix::matrix<DT_FP,char>* rot_matrix = _matrices_rotation[thread_number*MAX_DEPTH+depth];
						glLib::Math::Matrix::set_rotation_matrix_rad(rot_matrix,tangent,distribution_angular_offset);
						int num_valid_samples = 0;
						DT_FP avg_dist = 0.0;
						for (int i=0;i<trace_info.hit_material->samples_refraction;++i) {
							DT_FP random_cos_vec[3], perturbed_refraction[3];

							glLib::Math::Random::get_hemisphere_cosine_n( random_cos_vec, trace_info.hit_material->exponent_refractive, (DT_FP)(1.0) );
							for (char j=0;j<3;++j) _matrix_point->values[j][0] = random_cos_vec[j];

							glLib::Math::Matrix::multiply(
								rot_matrix,
								_matrix_point,
								_matrix_result
							);
							for (char j=0;j<3;++j) perturbed_refraction[j] = _matrix_result->get_at(0,j);

							//if (glLibMathVec::dot(perturbed_refraction,n,3)>0.0) {
							//	glLibMathVec::negate(perturbed_refraction,3);
							//	glLibMathVec::reflect(perturbed_refraction,n,perturbed_refraction,3);
							//}

							Ray refr_ray(offset_position,perturbed_refraction);
							RayInfo refract_result = trace_ray(refr_ray,refraction,depth+1,&material_stack_temp, thread_number);
							if (!refract_result.error) {
								++num_valid_samples;
								avg_dist += refract_result.hit_distance;
							}
							#ifdef ALLOW_BACKGROUND
							if (!refract_result.hit) {
								glLib::Math::Vec::add3(refraction,_renderer->_background_color,refraction);
							}
							#endif
						}
						avg_dist /= num_valid_samples;
						glLib::Math::Vec::scale3( (DT_COL)(1.0)/((DT_FP)(5.0)*avg_dist), refraction );
						glLib::Math::Vec::scale3( (DT_COL)(1.0)/(DT_COL)(num_valid_samples), refraction );
					}*/
					Ray ray_refraction(offset_position,vec_refraction);
					#ifdef DEBUG_CLICK
					if (_renderer->render_debug) {
						DEBUG_CLICK_INDENT printf("Tracing refraction ray:\n");
					}
					#endif
					RayInfo result_refract = trace_ray(ray_refraction,color_sample_refraction,depth+1,&material_stack_refract_ray);
					#ifdef DEBUG_CLICK
					if (_renderer->render_debug) {
						DEBUG_CLICK_INDENT printf("Done tracing refraction ray\n");
						DEBUG_CLICK_INDENT printf("-Refracted ray color: "); print_color(color_sample_refraction);
						//if (!result_refract.hit) printf(" + BG color");
						printf("\n");
					}
					#endif
					#ifdef ALLOW_BACKGROUND
					if (!result_refract.hit) {
						glLib::Math::Vec::add3(color_sample_refraction, _renderer->_background_color, color_sample_refraction);
					}
					#endif

					glLib::Math::Vec::multiply3(coefficient_transmission, color_sample_refraction,color_sample_refraction);
					glLib::Math::Vec::add3(color_sample_refraction, color,color);
					#ifdef DEBUG_CLICK
					if (_renderer->render_debug) {
						//DEBUG_CLICK_INDENT printf("-theta1=%3.2f; theta2=%3.2f\n",glLib::Math::Misc::degrees(theta_1),glLib::Math::Misc::degrees(theta_3));
						//DEBUG_CLICK_INDENT printf("-Weighting refraction sample by "); glLib::Math::Vec::print3(coefficient_transmission);
					}
					#endif
				}
			}
			#endif
			if (must_calculate_fresnel) {
				//Either total internal reflection, or material is not refractive
				//Set the Fresnel Coefficients, attenuating by thin films, if necessary (and allowed).
				#ifdef ALLOW_FILMS
				if (trace_info.hit_material->film) {
					thin_film_reflect( n1, theta_1, film_index,film_thickness, position_intersection, coefficient_reflection );
					//TODO:
					#ifdef OUTPUT_FILM_PATH_LENGTH
					add_color_heatmap(color, path_length/OUTPUT_FILM_PATH_LENGTH_SCALE); return trace_info;
					#endif
				} else {
				#endif
					glLib::Math::Vec::fillones3(coefficient_reflection);
				#ifdef ALLOW_FILMS
				}
				#endif
			}
			#ifdef ALLOW_REFLECTIONS
			if (trace_info.hit_material->reflect) {
				DT_COL color_sample_reflection[3] = {0.0,0.0,0.0};

				glLib::Math::Vec::reflect3(vec_source,vec_normal,vec_reflection);

				#ifdef OUTPUT_REFLECTION
				if (depth==OUTPUT_REFLECTION_DEPTH) {
					map_vec_to_rgb(r,color);
					return trace_info;
				}
				#endif

				DT_FP offset_position[3];
				glLib::Math::Vec::scaleadd3(EPSILON,vec_normal, position_intersection, offset_position);

				/*if (trace_info.hit_material->samples_reflection>1) {
					//TODO: define glLibMathVec::y_R3, x_R3, etc.
					DT_FP world_y[3] = {0.0,1.0,0.0};
					DT_FP tangent[3]; glLib::Math::Vec::cross(world_y,r,tangent); glLib::Math::Vec::normalize3(tangent);
					DT_FP distribution_angular_offset = acos(glLib::Math::Vec::dot3(world_y,r));

					glLib::Math::Matrix::matrix<DT_FP,char>* rot_matrix = _matrices_rotation[thread_number*MAX_DEPTH+depth];
					glLib::Math::Matrix::set_rotation_matrix_rad(rot_matrix,tangent,distribution_angular_offset);
					for (int i=0;i<trace_info.hit_material->samples_reflection;++i) {
						DT_FP random_cos_vec[3], perturbed_reflection[3];

						glLib::Math::Random::get_hemisphere_cosine_n(random_cos_vec,trace_info.hit_material->shininess,(DT_FP)(1.0));
						for (char j=0;j<3;++j) _matrix_point->values[j][0] = random_cos_vec[j];

						glLib::Math::Matrix::multiply(
							rot_matrix,
							_matrix_point,
							_matrix_result
						);
						for (char j=0;j<3;++j) perturbed_reflection[j] = _matrix_result->get_at(0,j);

						if (glLib::Math::Vec::dot3(perturbed_reflection,n)<0.0) {
							glLib::Math::Vec::negate3(perturbed_reflection);
							glLib::Math::Vec::reflect3(perturbed_reflection,n,perturbed_reflection);
						}

						deque<Material*> material_stack_temp = *material_stack;
						Ray refl_ray(offset_position,perturbed_reflection);
						RayInfo reflect_result = trace_ray(refl_ray,reflection,depth+1,&material_stack_temp, thread_number);
						#ifdef ALLOW_BACKGROUND
						if (!reflect_result.hit) {
							glLib::Math::Vec::add3(reflection,_renderer->_background_color,reflection);
						}
						#endif
					}

					glLib::Math::Vec::scale3((DT_COL)(1.0)/(DT_COL)(trace_info.hit_material->samples_reflection),reflection);
				}*/
				Ray ray_reflection(offset_position,vec_reflection);
				#ifdef DEBUG_CLICK
				if (_renderer->render_debug) {
					DEBUG_CLICK_INDENT printf("Tracing reflection ray:\n");
				}
				#endif
				RayInfo result_reflect = trace_ray(ray_reflection,color_sample_reflection,depth+1,&material_stack_reflect_ray);
				#ifdef DEBUG_CLICK
				if (_renderer->render_debug) {
					DEBUG_CLICK_INDENT printf("Done tracing reflection ray\n");
					DEBUG_CLICK_INDENT printf("-Reflected ray color: "); print_color(color_sample_reflection);
					//if (!result_reflect.hit) printf(" + BG color");
					printf("\n");
				}
				#endif
				#ifdef ALLOW_BACKGROUND
				if (!result_reflect.hit) {
					glLib::Math::Vec::add3(color_sample_reflection, _renderer->_background_color, color_sample_reflection);
				}
				#endif

				glLib::Math::Vec::multiply3(coefficient_reflection, color_sample_reflection,color_sample_reflection);
				glLib::Math::Vec::add3(color_sample_reflection, color,color);
				#ifdef DEBUG_CLICK
				if (_renderer->render_debug) {
					//DEBUG_CLICK_INDENT printf("-Weighting reflection sample by "); glLib::Math::Vec::print3(coefficient_reflection);
				}
				#endif
			}
			#endif

			#undef film_index
			#undef film_thickness

			#ifdef OUTPUT_FRESNEL
			color[0] = coefficient_reflection[0];
			color[1] = 0.0;
			color[2] = coefficient_transmission[1];
			return trace_info;
			#endif
		}
	}

	//Attenuate
	color[0] *= exp(-material_current->absorption[0]*trace_info.hit_distance);
	color[1] *= exp(-material_current->absorption[1]*trace_info.hit_distance);
	color[2] *= exp(-material_current->absorption[2]*trace_info.hit_distance);
	/*if (material_current->absorption[0]!=0.0) {
		color[0] += trace_info.hit_distance*0.1;
		color[1] += trace_info.hit_distance*0.1;
		color[2] += trace_info.hit_distance*0.1;
	}*/
	//for (int i=0;i<depth;++i){printf("  ");} printf("Traveling through material %s\n",material_current->name->c_str());

	#ifdef DEBUG_CLICK
	if (_renderer->render_debug) {
		//DEBUG_CLICK_INDENT printf("-Attenuating color by distance %e\n",trace_info.hit_distance);
	}
	#endif

	return trace_info;
}
#else
RayInfo TileRenderer::trace_shadow_ray(Ray ray) {}
RayInfo TileRenderer::trace_ray(Ray ray_incident, DT_COL color[3], int depth, deque<Material*>* material_stack) {
	RayInfo trace_info;
	unsigned int num_tests = trace_get_intersect(ray_incident,&trace_info);
	if (!trace_info.hit) return trace_info;

	Object* obj = trace_info.hit_object;

	//Get source vector direction
	DT_FP vec_source[3];
	glLib::Math::Vec::copy3(ray_incident.get_direction_negative(),vec_source);

	//Get intersection point
	DT_FP position_intersection[3];
	//scaleadd TODO
	glLib::Math::Vec::scale3(trace_info.hit_distance,ray_incident.get_direction(),position_intersection);
	glLib::Math::Vec::add3(position_intersection,ray_incident.origin,position_intersection);

	//Get normal
	DT_FP vec_normal[3];
	trace_info.hit_object->set_normal(position_intersection,vec_normal);
	if (glLib::Math::Vec::dot3(vec_source,vec_normal)<0.0) {
		//Change the normal so as to be facing toward the incident vector.
		glLib::Math::Vec::negate3(vec_normal);
	}

	//Offset intersection position
	glLib::Math::Vec::scaleadd3(EPSILON,vec_normal, position_intersection, position_intersection);

	DT_COL f[3]; glLib::Math::Vec::copy3(obj->material->diffuse,f);

	if (depth>5) {
		double p = max(f[0],max(f[1],f[2]));
		if (glLib::Math::Random::get_uniform<double>()<p) glLib::Math::Vec::scale3(1.0/p,f);
		else {
			glLib::Math::Vec::add3(obj->material->emission,color, color);
			return trace_info;
		}
	}

	DT_COL sample[3] = {0.0,0.0,0.0}; DT_FP ray_direction[3];

	//Ideal diffuse reflection
	glLib::Math::Random::get_hemisphere_cosine(ray_direction,vec_normal);
	
	//Ideal specular reflection
	//glLib::Math::Vec::reflect3(vec_source,vec_normal,ray_direction);
	
/*
	// Ideal dielectric REFRACTION
	Ray reflRay(x, r.d-n*2*n.dot(r.d));
	bool into = n.dot(nl)>0; // Ray from outside going in?
	double nc = 1;
	double nt = 1.5;
	double nnt=into?nc/nt:nt/nc;
	double ddn = r.d.dot(nl);
	double cos2t;
	if ((cos2t=1-nnt*nnt*(1-ddn*ddn))<0) { // Total internal reflection
		return obj.e + f.mult(radiance(reflRay,depth));
	}

	Vec tdir = (r.d*nnt - n*((into?1:-1)*(ddn*nnt+sqrt(cos2t)))).norm();
	double a = nt - nc;
	double b = nt + nc;
	double R0 = (a*a) / (b*b);
	double c = 1.0 - (into ? -ddn : tdir.dot(n));
	double Re = R0+(1-R0)*c*c*c*c*c;
	double Tr = 1.0 - Re;
	double P = 0.25 + 0.5*Re;
	double RP = Re / P;
	double TP = Tr / (1.0-P);

	// Russian roulette
	return obj.e + f.mult(
		depth>2 ?
			(runiform()<P ?
				radiance(reflRay,depth)*RP
				:
				radiance(Ray(x,tdir),depth)*TP
			)
			:
			radiance(reflRay,depth) * Re + radiance(Ray(x,tdir),depth)*Tr
	);*/

	Ray ray(position_intersection,ray_direction);
	trace_ray(ray,sample,depth+1,material_stack);

	glLib::Math::Vec::multiply3(f,sample,sample);
	glLib::Math::Vec::add3(sample,color, color);
	glLib::Math::Vec::add3(obj->material->emission,color, color);

	return trace_info;
}
#endif
Uint32 TileRenderer::generate_trace_rays(int pix_x, int pix_y, DT_COL color_sample[3], deque<Material*>* material_stack) {
	glLib::Math::Vec::fillzeros3(color_sample);
	for (int i=0;i<_renderer->_numof_samples_pixel;++i) {
		material_stack->clear();
		material_stack->push_front(material_air);

		DT_FP x_part = ((DT_FP)(pix_x)+_renderer->_poisson_samples_pixel[i]->x) / (DT_FP)_renderer->_img_w * -2.0 + 1.0;
		DT_FP y_part = ((DT_FP)(pix_y)+_renderer->_poisson_samples_pixel[i]->y) / (DT_FP)_renderer->_img_h * -2.0 + 1.0;
		y_part /= _renderer->img_aspect;

		DT_FP ray_position[3] = {
			x_part*_renderer->vec_eye_x[0] + y_part*_renderer->vec_eye_y[0] + _renderer->vec_img_center[0],
			x_part*_renderer->vec_eye_x[1] + y_part*_renderer->vec_eye_y[1] + _renderer->vec_img_center[1],
			x_part*_renderer->vec_eye_x[2] + y_part*_renderer->vec_eye_y[2] + _renderer->vec_img_center[2]
		};

		//Calculate focal point of this sample in the scene
		// 1/s1 + 1/s2 = 1/f
		//s2 is _scene->camera->image_plane
		//s1 is the camera's z distance from the lens to the focal point
		DT_FP central_ray[3], central_ray_prime[3]; //Rays through the lens' center
		glLib::Math::Vec::subtract3(_renderer->vec_camera_position,ray_position,central_ray);
		static DT_FP s2 = _renderer->_scene->camera->image_plane;
		static DT_FP s1 = _renderer->_scene->camera->recip_focal_length - 1.0/s2;
		//central_ray_prime is parallel to central_ray, but goes exactly to the focal point
		// ||central_ray||/s1 = ||central_ray_prime||/(s1+s2)
		glLib::Math::Vec::scale3((s1+s2)*glLib::Math::Vec::length3(central_ray)/s1,central_ray,central_ray_prime);
		DT_FP focal_point[3]; glLib::Math::Vec::add3(ray_position,central_ray_prime,focal_point);

		int distr_index = _renderer->_numof_samples_aperture*(rand()%NUM_RANDOM_DISTR);
		for (int j=0;j<_renderer->_numof_samples_aperture;++j) {
			DT_FP ray_direction[3]; //ray_direction is also used as a temp
			glLib::Math::Vec::scale3(_renderer->_poisson_samples_aperture[distr_index+j]->x,_renderer->vec_eye_x,ray_direction); //ray_direction is used as temp
			glLib::Math::Vec::add3(_renderer->vec_camera_position,ray_direction,ray_position);
			glLib::Math::Vec::scale3(_renderer->_poisson_samples_aperture[distr_index+j]->y,_renderer->vec_eye_y,ray_direction); //ray_direction is used as temp
			glLib::Math::Vec::add3(ray_direction,ray_position,ray_position);

			//Ray starts on the lens, heads toward focal point
			glLib::Math::Vec::subtract3(focal_point,ray_position,ray_direction); //ray_direction gets set for real
			glLib::Math::Vec::normalize3(ray_direction);

			#ifdef DEBUG_CLICK
			if (_renderer->render_debug) printf("Begin tracing ray:\n");
			#endif

			Ray ray(ray_position,ray_direction);
			#ifndef ALLOW_BACKGROUND
			trace_ray(ray,color,0,material_stack);
			#else
			if (!trace_ray(ray,color_sample,0,material_stack).hit) {
				glLib::Math::Vec::add3(color_sample,_renderer->_background_color,color_sample);
			}
			#endif

			#ifdef DEBUG_CLICK
			if (_renderer->render_debug) printf("Traced ray!\n");
			#endif
		}
	}

	#ifdef DEBUG_CLICK
	if (_renderer->render_debug) {
		printf("Pixel color: "); print_color(color_sample); printf("\n");
	}
	#endif

	//[0.0,n*1.0] => [0.0,255.0]
	glLib::Math::Vec::scale3(
		(DT_COL)(255.0) / (_renderer->_numof_samples_pixel*_renderer->_numof_samples_aperture),
		color_sample
	);
	glLib::Math::Vec::clamp3(color_sample,(DT_COL)(0.0),(DT_COL)(255.0));

	Uint32 final_color = SDL_MapRGB(
		_renderer->_pixel_format,
		glLib::Math::Misc::rndchar(color_sample[0]),
		glLib::Math::Misc::rndchar(color_sample[1]),
		glLib::Math::Misc::rndchar(color_sample[2])
	);

	return final_color;
}
int TileRenderer::render_tile(void* data) {
	deque<Material*> material_stack = deque<Material*>();

	#ifdef SAVE_INTERMEDIARY
	string tile_name;
	SDL_Surface* current_tile = SDL_CreateRGBSurface(SDL_SWSURFACE,TILE_SIZE,TILE_SIZE,24,0,0,0,0);
	#endif
	while (true) {
		Tile* to_render = _renderer->get_next_tile();
		if (to_render==NULL) break;

		#ifdef ALLOW_RENDER
		for (int pixel_x=to_render->x;pixel_x<to_render->x+to_render->w;pixel_x+=PIXEL_STEP) {
			for (int pixel_y=to_render->y;pixel_y<to_render->y+to_render->h;pixel_y+=PIXEL_STEP) {
				DT_COL color_sample[3];
				
				Uint32 final_color = generate_trace_rays(pixel_x,pixel_y,color_sample,&material_stack);
				
				set_pixel_rgb24(pixel_x,_renderer->_img_h-pixel_y-1,final_color,_renderer->_image);
				#ifdef SAVE_INTERMEDIARY
				set_pixel_rgb24(pixel_x-to_render.x,TILE_SIZE-(pixel_y-to_render.y)-1,final_color,current_tile);
				#endif
			}
		}
		#else
		SDL_Delay(RENDER_TILE_DELAY);
		#endif

		#ifdef SAVE_INTERMEDIARY
		save(tile_name.c_str(),current_tile);
		#endif

		_renderer->complete_tile(to_render);
	}

	#ifdef SAVE_INTERMEDIARY
	SDL_FreeSurface(current_tile);
	#endif

	return 0;
}