Index: GTP/trunk/Lib/Vis/Preprocessing/src/GvsPreprocessor.h
===================================================================
--- GTP/trunk/Lib/Vis/Preprocessing/src/GvsPreprocessor.h	(revision 2738)
+++ GTP/trunk/Lib/Vis/Preprocessing/src/GvsPreprocessor.h	(revision 2739)
@@ -121,17 +121,5 @@
 
 protected:
-#if 0
-	struct PendingRay
-	{
-		PendingRay(VssRay *ray, const bool d)
-			: mRay(ray), mFoundDiscontinuity(d) 
-		{}
-
-		VssRay *mRay;
-		bool mFoundDiscontinuity;
-	};
-
-	typedef stack<PendingRay> PendingQueue;
-#endif
+
 	typedef stack<VssRay *> RayQueue;
 
@@ -140,13 +128,10 @@
 		by the queue and continues as long it finds new visible geometry 
 		(i.e., the queue is	not empty).
-
 		@returns the number of samples cast.
 	*/
 	int ProcessQueue();
-
 	/** Generates the rays starting the adaptive visibility sampling process.
 	*/
 	int CastInitialSamples(int numSamples);
-
 	/** Uses the information gained from the ray for doing adaptive border sampling.
 		This function tries to find the border of the triangle found visible by the
@@ -163,4 +148,5 @@
 	*/	
 	int AdaptiveBorderSampling(const VssRay &currentRay);
+
 	int AdaptiveBorderSamplingOpt(const VssRay &currentRay);
 	/** The reverse sampling step. It is started once the cast
@@ -179,4 +165,10 @@
 						   const Triangle3 &hitTriangle,
 						   const VssRay &oldRay);
+	/** Returns true if we sampled a closer triangle than with the previous ray. 
+		Does reverse sampling if gap found.
+	*/
+	int CheckDiscontinuity2(const VssRay &currentRay, 
+		                    const Triangle3 &hitTriangle,
+							const VssRay &oldRay);
 	/** Adds new samples to the ray queue and classifies them
 		with respect to the previous ray.
@@ -300,4 +292,17 @@
 
 	void PrepareProbablyVisibleSampling();
+	
+	void CreateRandomizedReverseRays(const Vector3 &origin, 
+		                             const Vector3 &termination, 
+									 const Triangle3 &triangle,                   
+									 SimpleRayContainer &rays, 
+									 int number);
+
+	bool CreateRandomizedReverseRay(const Vector3 &origin, 
+		                            const Vector3 &termination,
+									const Triangle3 &triangle,
+									SimpleRay &ray);
+
+	inline bool IntersectViewCell(Vector3 &origin, Vector3 &dir) const;
 
 
Index: GTP/trunk/Lib/Vis/Preprocessing/src/GvsPreprocessor.cpp
===================================================================
--- GTP/trunk/Lib/Vis/Preprocessing/src/GvsPreprocessor.cpp	(revision 2738)
+++ GTP/trunk/Lib/Vis/Preprocessing/src/GvsPreprocessor.cpp	(revision 2739)
@@ -176,7 +176,6 @@
 	// discrepancy between predicted and actual hit point large => this is likely to be a discontinuity
 
-	// matt: problem: this does only find CLOSER discontinuities, but not discontinuities which are FARTHER away
+	// matt: problem: this only finds CLOSER discontinuities, but not discontinuities which are FARTHER away
 	if ((predictedLen - len) < mThreshold)
-	//if ((fabs(predictedLen - len) < mThreshold)
 	// q: rather use relative distance?
 	//if ((predictedLen / len) < mThreshold)
@@ -208,4 +207,36 @@
 			CastRayBundle16(simpleRay, reverseRays, scale);
 	}
+
+	mGvsStats.mReverseSamples += (int)reverseRays.size();
+
+	EnqueueRays(reverseRays);
+
+	return (int)reverseRays.size();
+}
+
+
+int GvsPreprocessor::CheckDiscontinuity2(const VssRay &currentRay,
+										 const Triangle3 &hitTriangle,
+										 const VssRay &oldRay)
+{
+	// the predicted hitpoint is the point where the ray would have intersected the hit triangle
+	Vector3 predictedHit = CalcPredictedHitPoint(currentRay, hitTriangle, oldRay);
+
+	float predictedLen = Magnitude(predictedHit - currentRay.mOrigin);
+	float len = Magnitude(currentRay.mTermination - currentRay.mOrigin);
+	
+	// discrepancy between predicted and actual hit point large 
+	// => this is likely to be a discontinuity
+	if (fabs(predictedLen - len) < mThreshold)
+		return 0;
+
+	static SimpleRayContainer simpleRays;
+	simpleRays.clear();
+
+	CreateRandomizedReverseRays(currentRay.mOrigin, predictedHit, hitTriangle, simpleRays, 16);
+
+	static VssRayContainer reverseRays;
+	reverseRays.clear();
+	CastRays(simpleRays, reverseRays, false, false);
 
 	mGvsStats.mReverseSamples += (int)reverseRays.size();
@@ -487,4 +518,38 @@
 	}
 
+	//////////
+	//-- fill up simple rays with random rays so we can cast 16 rays at a time
+
+	const int numRandomRays = 16 - (int)simpleRays.size();
+	GenerateRays(numRandomRays, *mDistribution, simpleRays, sInvalidSamples);
+
+	generationTimer.Exit();
+
+
+	/////////////////////
+	//-- cast rays to vertices and determine visibility
+	
+	static VssRayContainer vssRays;
+	vssRays.clear();
+
+	// for the original implementation we don't cast double rays
+	const bool castDoubleRays = !mPerViewCell;
+	// cannot prune invalid rays because we have to compare adjacent rays
+	const bool pruneInvalidRays = false;
+
+	//gvsTimer.Entry();
+	castTimer.Entry();
+
+	CastRays(simpleRays, vssRays, castDoubleRays, pruneInvalidRays);
+	
+	castTimer.Exit();
+	//gvsTimer.Exit();
+
+	const int numBorderSamples = (int)vssRays.size() - numRandomRays;
+	int castRays = (int)simpleRays.size();
+
+	// handle cast rays
+	EnqueueRays(vssRays);
+
 	if (0)
 	{
@@ -493,60 +558,7 @@
 		dummy.enlargedTriangle = new Polygon3(enlargedTriangle);
 		dummy.originalTriangle = hitTriangle;
+	
 		vizContainer.push_back(dummy);
-	}
-
-	
-	//////////
-	//-- fill up simple rays with random rays so we can cast 16 rays at a time
-
-	const int numRandomRays = 16 - (int)simpleRays.size();
-
-	if (!mUseDeterministicGvs)
-	{
-		if (0)
-		{
-			SimpleRay mainRay = SimpleRay(currentRay.GetOrigin(), 
-				                          currentRay.GetNormalizedDir(), 
-										  SamplingStrategy::GVS, 
-										  1.0f);
-			
-			GenerateJitteredRays(simpleRays, mainRay, numRandomRays, 0, 0.01f);
-		}
-		else
-			GenerateImportanceSamples(currentRay, hitTriangle, numRandomRays, simpleRays);
-	}
-	else
-		GenerateRays(numRandomRays, *mDistribution, simpleRays, sInvalidSamples);
-
-	generationTimer.Exit();
-
-
-	/////////////////////
-	//-- cast rays to vertices and determine visibility
-	
-	static VssRayContainer vssRays;
-	vssRays.clear();
-
-	// for the original implementation we don't cast double rays
-	const bool castDoubleRays = !mPerViewCell;
-	// cannot prune invalid rays because we have to compare adjacent rays
-	const bool pruneInvalidRays = false;
-
-	//gvsTimer.Entry();
-	castTimer.Entry();
-
-	CastRays(simpleRays, vssRays, castDoubleRays, pruneInvalidRays);
-	
-	castTimer.Exit();
-	//gvsTimer.Exit();
-
-	const int numBorderSamples = (int)vssRays.size() - numRandomRays;
-	int castRays = (int)simpleRays.size();
-
-	// handle cast rays
-	EnqueueRays(vssRays);
-
-	if (0)
-	{
+	
 		// set flags
 		VssRayContainer::const_iterator rit, rit_end = vssRays.end();
@@ -556,14 +568,9 @@
 	}
 
-	for (size_t i = 0; i < vssRays.size(); ++ i)
-	{
-		// check for discontinuity and cast reverse rays if necessary
-		castRays += CheckDiscontinuity(*vssRays[i], hitTriangle, currentRay);
-	}
-	
-
     // recursivly subdivide each edge
 	for (int i = 0; i < numBorderSamples; ++ i)
 	{
+		castRays += CheckDiscontinuity(*vssRays[i], hitTriangle, currentRay);
+
 		castRays += SubdivideEdge(hitTriangle,
 								  enlargedTriangle[i], 
@@ -638,21 +645,6 @@
 	const int numRandomRays = 16 - ((int)simpleRays.size() % 16);
 
-		if (!mUseDeterministicGvs)
-	{
-		if (0)
-		{
-			SimpleRay mainRay = SimpleRay(currentRay.GetOrigin(), 
-				                          currentRay.GetNormalizedDir(), 
-										  SamplingStrategy::GVS, 
-										  1.0f);
-			
-			GenerateJitteredRays(simpleRays, mainRay, numRandomRays, 0, 0.01f);
-		}
-		else
-			GenerateImportanceSamples(currentRay, hitTriangle, numRandomRays, simpleRays);
-	}
-	else
-		GenerateRays(numRandomRays, *mDistribution, simpleRays, sInvalidSamples);
-
+	GenerateImportanceSamples(currentRay, hitTriangle, numRandomRays, simpleRays);
+	
 	generationTimer.Exit();
 
@@ -689,5 +681,5 @@
 	{
 		// check for discontinuity and cast reverse rays if necessary
-		castRays += CheckDiscontinuity(*vssRays[i], hitTriangle, currentRay);
+		castRays += CheckDiscontinuity2(*vssRays[i], hitTriangle, currentRay);
 	}
 
@@ -794,26 +786,35 @@
 
 	//////////////
-	//-- for per view cell sampling, we must check for intersection
+	//-- for per view cell sampling, we must intersect
 	//-- with the current view cell
 
     if (mPerViewCell)
 	{
-		// send ray to view cell
-		static Ray ray;
-
-		ray.Clear();
-		ray.Init(newOrigin, -newDir, Ray::LOCAL_RAY);
-		
-		// check if ray intersects view cell
-		if (!mCurrentViewCell->CastRay(ray))
+		if (!IntersectViewCell(newOrigin, newDir))
 			return false;
-
-		Ray::Intersection &hit = ray.intersections[0];
-	
-		// the ray starts from the view cell
-		newOrigin = ray.Extrap(hit.mT);
 	}
 
 	reverseRay = SimpleRay(newOrigin, newDir, SamplingStrategy::GVS, 1.0f);
+
+	return true;
+}
+
+
+bool GvsPreprocessor::IntersectViewCell(Vector3 &origin, Vector3 &dir) const
+{
+	// send ray to view cell
+	static Ray testRay;
+
+	testRay.Clear();
+	testRay.Init(origin, -dir, Ray::LOCAL_RAY);
+
+	// check if ray intersects view cell
+	if (!mCurrentViewCell->CastRay(testRay))
+		return false;
+
+	Ray::Intersection &hit = testRay.intersections[0];
+
+	// the ray starts from the view cell
+	origin = testRay.Extrap(hit.mT);
 
 	return true;
@@ -908,8 +909,9 @@
 		mRayQueue.pop();
 		
-		if (0)
+		if (mUseDeterministicGvs)
 			castSamples += AdaptiveBorderSampling(*ray);
 		else
 			castSamples += AdaptiveBorderSamplingOpt(*ray);
+
 		// note: don't have to delete because handled by ray pool
 	}
@@ -2033,3 +2035,70 @@
 
 
-}
+void GvsPreprocessor::CreateRandomizedReverseRays(const Vector3 &origin, 
+												  const Vector3 &termination, 
+												  const Triangle3 &triangle,                   
+												  SimpleRayContainer &rays, 
+												  int number)
+{
+	size_t currentNumber = rays.size();
+
+	while ((rays.size() - currentNumber) < number)
+	{
+		SimpleRay ray;
+		if (CreateRandomizedReverseRay(origin, termination, triangle, ray))
+			rays.push_back(ray);
+	}
+}
+
+
+bool GvsPreprocessor::CreateRandomizedReverseRay(const Vector3 &origin, 
+												 const Vector3 &termination, 
+												 const Triangle3 &triangle,                   
+												 SimpleRay &ray)
+{
+	Vector3 nd = triangle.GetNormal();
+  
+	// Compute right handed coordinate system from direction
+	Vector3 U, V;
+	nd.RightHandedBase(U, V);
+	
+	float rr[2];
+
+	rr[0] = RandomValue(0, 1);
+	rr[1] = RandomValue(0, 1);
+
+	Vector2 vr2(rr[0], rr[1]);
+	const float sigma = triangle.GetBoundingBox().Radius() * mRadiusOfInfluence;
+	Vector2 gaussvec2;
+
+	GaussianOn2D(vr2, sigma, // input
+		         gaussvec2); // output
+	
+	const Vector3 shift = gaussvec2.xx * U + gaussvec2.yy * V;
+
+	Vector3 newOrigin = origin;
+
+	//cout << "t: " << termination;
+	newOrigin += shift;
+	//cout << " new t: " << termination << endl;
+	Vector3 direction = termination - newOrigin;
+
+	const float len = Magnitude(direction);
+	
+	if (len < Limits::Small)
+		return false;
+  
+	direction /= len;
+
+	if (!IntersectViewCell(newOrigin, direction))
+		return false;
+
+	// $$ jb the pdf is yet not correct for all sampling methods!
+	const float pdf = 1.0f;
+	ray = SimpleRay(newOrigin, direction, SamplingStrategy::GVS, pdf);
+	
+	return true;
+}
+
+
+}
