 	// set texture parameters
 	tex.normalized = true;                      // access with normalized texture coordinates
 	tex.filterMode = cudaFilterModeLinear;      // linear interpolation
-	tex.addressMode[0] = cudaAddressModeClamp; // we want 0s outside for beautiful isosurfaces at the borders
-	tex.addressMode[1] = cudaAddressModeClamp;
+	tex.addressMode[0] = cudaAddressModeMirror;
+	tex.addressMode[1] = cudaAddressModeMirror;
 
 	// bind array to 3D texture
 	checkCudaErrors(cudaBindTextureToArray(tex, d_volumeArray, channelDesc));

#include <thrust/extrema.h>
#include <thrust/device_ptr.h>

typedef unsigned int  uint;
typedef unsigned char uchar;


		samplePos.x = (pos.x - boxMin.x) / diff.x;
		samplePos.y = (pos.y - boxMin.y) / diff.y;
		samplePos.z = (pos.z - boxMin.z) / diff.z;

		float curDist = 0.0f;
#ifdef CYLINDER
		curDist = sqrt((0.5f - samplePos.x) * (0.5f - samplePos.x) + (0.5f - samplePos.y) * (0.5f - samplePos.y));
		if (curDist > 0.5) {
			float2 middir = make_float2(0.5f, 0.5f) - make_float2(samplePos.x, samplePos.y);
			middir = normalize(middir);
			samplePos.x = samplePos.x + (curDist - 0.5f) * middir.x;
			samplePos.y = samplePos.y + (curDist - 0.5f) * middir.y;
		}
#endif CYLINDER

		float sample = (tex3D(tex, samplePos.x, samplePos.y, samplePos.z) - minVal) / (maxVal - minVal);

		if (pointBehindPlane(samplePos, plane) || curDist > 0.5f) {

			for (int isoIndex = 0; isoIndex < d_numIsoVals; isoIndex++) {
				isoDiffs[isoIndex] = sample - isoVals[isoIndex];
