Merge branch 'texture'

This replaces the deprecated CUDA texture reference API by texture objects.

10 files changed, 275 insertions, 291 deletions

diff --git a/cuda/2d/fan_bp.cu b/cuda/2d/fan_bp.cu
index 2068d03..bcda464 100644
--- a/cuda/2d/fan_bp.cu
+++ b/cuda/2d/fan_bp.cu
@@ -57,25 +57,8 @@ struct DevFanParams {
 
 __constant__ DevFanParams gC_C[g_MaxAngles];
 
-
-static bool bindProjDataTexture(float* data, unsigned int pitch, unsigned int width, unsigned int height, cudaTextureAddressMode mode = cudaAddressModeBorder)
-{
-	cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<float>();
-
-	gT_FanProjTexture.addressMode[0] = mode;
-	gT_FanProjTexture.addressMode[1] = mode;
-	gT_FanProjTexture.filterMode = cudaFilterModeLinear;
-	gT_FanProjTexture.normalized = false;
-
-	cudaBindTexture2D(0, gT_FanProjTexture, (const void*)data, channelDesc, width, height, sizeof(float)*pitch);
-
-	// TODO: error value?
-
-	return true;
-}
-
 template<bool FBPWEIGHT>
-__global__ void devFanBP(float* D_volData, unsigned int volPitch, unsigned int startAngle, const SDimensions dims, float fOutputScale)
+__global__ void devFanBP(float* D_volData, unsigned int volPitch, cudaTextureObject_t tex, unsigned int startAngle, const SDimensions dims, float fOutputScale)
 {
 	const int relX = threadIdx.x;
 	const int relY = threadIdx.y;
@@ -115,7 +98,7 @@ __global__ void devFanBP(float* D_volData, unsigned int volPitch, unsigned int s
 
 		const float fr = __fdividef(1.0f, fDen);
 		const float fT = fNum * fr;
-		fVal += tex2D(gT_FanProjTexture, fT, fA) * (FBPWEIGHT ? fr * fr : fr);
+		fVal += tex2D<float>(tex, fT, fA) * (FBPWEIGHT ? fr * fr : fr);
 		fA += 1.0f;
 	}
 
@@ -123,7 +106,7 @@ __global__ void devFanBP(float* D_volData, unsigned int volPitch, unsigned int s
 }
 
 // supersampling version
-__global__ void devFanBP_SS(float* D_volData, unsigned int volPitch, unsigned int startAngle, const SDimensions dims, float fOutputScale)
+__global__ void devFanBP_SS(float* D_volData, unsigned int volPitch, cudaTextureObject_t tex, unsigned int startAngle, const SDimensions dims, float fOutputScale)
 {
 	const int relX = threadIdx.x;
 	const int relY = threadIdx.y;
@@ -169,7 +152,7 @@ __global__ void devFanBP_SS(float* D_volData, unsigned int volPitch, unsigned in
 				const float fr = __fdividef(1.0f, fDen);
 
 				const float fT = fNum * fr;
-				fVal += tex2D(gT_FanProjTexture, fT, fA) * fr;
+				fVal += tex2D<float>(tex, fT, fA) * fr;
 				fY -= fSubStep;
 			}
 			fX += fSubStep;
@@ -184,7 +167,7 @@ __global__ void devFanBP_SS(float* D_volData, unsigned int volPitch, unsigned in
 // BP specifically for SART.
 // It includes (free) weighting with voxel weight.
 // It assumes the proj texture is set up _without_ padding, unlike regular BP.
-__global__ void devFanBP_SART(float* D_volData, unsigned int volPitch, const SDimensions dims, float fOutputScale)
+__global__ void devFanBP_SART(float* D_volData, unsigned int volPitch, cudaTextureObject_t tex, const SDimensions dims, float fOutputScale)
 {
 	const int relX = threadIdx.x;
 	const int relY = threadIdx.y;
@@ -213,7 +196,7 @@ __global__ void devFanBP_SART(float* D_volData, unsigned int volPitch, const SDi
 	const float fr = __fdividef(1.0f, fDen);
 	const float fT = fNum * fr;
 	// NB: The scale constant in devBP is cancelled out by the SART weighting
-	const float fVal = tex2D(gT_FanProjTexture, fT, 0.5f);
+	const float fVal = tex2D<float>(tex, fT, 0.5f);
 
 	volData[Y*volPitch+X] += fVal * fOutputScale;
 }
@@ -303,11 +286,15 @@ bool FanBP_internal(float* D_volumeData, unsigned int volumePitch,
 {
 	assert(dims.iProjAngles <= g_MaxAngles);
 
-	bindProjDataTexture(D_projData, projPitch, dims.iProjDets, dims.iProjAngles);
+	cudaTextureObject_t D_texObj;
+	if (!createTextureObjectPitch2D(D_projData, D_texObj, projPitch, dims.iProjDets, dims.iProjAngles))
+		return false;
 
 	bool ok = transferConstants(angles, dims.iProjAngles, false);
-	if (!ok)
+	if (!ok) {
+		cudaDestroyTextureObject(D_texObj);
 		return false;
+	}
 
 	dim3 dimBlock(g_blockSlices, g_blockSliceSize);
 	dim3 dimGrid((dims.iVolWidth+g_blockSlices-1)/g_blockSlices,
@@ -318,15 +305,17 @@ bool FanBP_internal(float* D_volumeData, unsigned int volumePitch,
 
 	for (unsigned int i = 0; i < dims.iProjAngles; i += g_anglesPerBlock) {
 		if (dims.iRaysPerPixelDim > 1)
-			devFanBP_SS<<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, i, dims, fOutputScale);
+			devFanBP_SS<<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, D_texObj, i, dims, fOutputScale);
 		else
-			devFanBP<false><<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, i, dims, fOutputScale);
+			devFanBP<false><<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, D_texObj, i, dims, fOutputScale);
 	}
 
 	ok = checkCuda(cudaStreamSynchronize(stream), "FanBP");
 
 	cudaStreamDestroy(stream);
 
+	cudaDestroyTextureObject(D_texObj);
+
 	return ok;
 }
 
@@ -337,11 +326,15 @@ bool FanBP_FBPWeighted_internal(float* D_volumeData, unsigned int volumePitch,
 {
 	assert(dims.iProjAngles <= g_MaxAngles);
 
-	bindProjDataTexture(D_projData, projPitch, dims.iProjDets, dims.iProjAngles);
+	cudaTextureObject_t D_texObj;
+	if (!createTextureObjectPitch2D(D_projData, D_texObj, projPitch, dims.iProjDets, dims.iProjAngles))
+		return false;
 
 	bool ok = transferConstants(angles, dims.iProjAngles, true);
-	if (!ok)
+	if (!ok) {
+		cudaDestroyTextureObject(D_texObj);
 		return false;
+	}
 
 	dim3 dimBlock(g_blockSlices, g_blockSliceSize);
 	dim3 dimGrid((dims.iVolWidth+g_blockSlices-1)/g_blockSlices,
@@ -351,13 +344,15 @@ bool FanBP_FBPWeighted_internal(float* D_volumeData, unsigned int volumePitch,
 	cudaStreamCreate(&stream);
 
 	for (unsigned int i = 0; i < dims.iProjAngles; i += g_anglesPerBlock) {
-		devFanBP<true><<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, i, dims, fOutputScale);
+		devFanBP<true><<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, D_texObj, i, dims, fOutputScale);
 	}
 
 	ok = checkCuda(cudaStreamSynchronize(stream), "FanBP_FBPWeighted");
 
 	cudaStreamDestroy(stream);
 
+	cudaDestroyTextureObject(D_texObj);
+
 	return ok;
 }
 
@@ -369,19 +364,27 @@ bool FanBP_SART(float* D_volumeData, unsigned int volumePitch,
                 float fOutputScale)
 {
 	// only one angle
-	bindProjDataTexture(D_projData, projPitch, dims.iProjDets, 1, cudaAddressModeClamp);
+	cudaTextureObject_t D_texObj;
+	if (!createTextureObjectPitch2D(D_projData, D_texObj, projPitch, dims.iProjDets, 1, cudaAddressModeClamp))
+		return false;
 
 	bool ok = transferConstants(angles + angle, 1, false);
-	if (!ok)
+	if (!ok) {
+		cudaDestroyTextureObject(D_texObj);
 		return false;
+	}
 
 	dim3 dimBlock(g_blockSlices, g_blockSliceSize);
 	dim3 dimGrid((dims.iVolWidth+g_blockSlices-1)/g_blockSlices,
 	             (dims.iVolHeight+g_blockSliceSize-1)/g_blockSliceSize);
 
-	devFanBP_SART<<<dimGrid, dimBlock>>>(D_volumeData, volumePitch, dims, fOutputScale);
+	devFanBP_SART<<<dimGrid, dimBlock>>>(D_volumeData, volumePitch, D_texObj, dims, fOutputScale);
+
+	ok = checkCuda(cudaThreadSynchronize(), "FanBP_SART");
 
-	return checkCuda(cudaThreadSynchronize(), "FanBP_SART");
+	cudaDestroyTextureObject(D_texObj);
+
+	return ok;
 }
 
 bool FanBP(float* D_volumeData, unsigned int volumePitch,
diff --git a/cuda/2d/fan_fp.cu b/cuda/2d/fan_fp.cu
index 342ca4c..5aa9674 100644
--- a/cuda/2d/fan_fp.cu
+++ b/cuda/2d/fan_fp.cu
@@ -33,12 +33,6 @@ along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
 #include <iostream>
 #include <list>
 
-
-typedef texture<float, 2, cudaReadModeElementType> texture2D;
-
-static texture2D gT_FanVolumeTexture;
-
-
 namespace astraCUDA {
 
 static const unsigned g_MaxAngles = 2560;
@@ -55,31 +49,9 @@ static const unsigned int g_anglesPerBlock = 16;
 static const unsigned int g_detBlockSize = 32;
 static const unsigned int g_blockSlices = 64;
 
-static bool bindVolumeDataTexture(float* data, cudaArray*& dataArray, unsigned int pitch, unsigned int width, unsigned int height)
-{
-	cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<float>();
-	dataArray = 0;
-	cudaMallocArray(&dataArray, &channelDesc, width, height);
-	cudaMemcpy2DToArray(dataArray, 0, 0, data, pitch*sizeof(float), width*sizeof(float), height, cudaMemcpyDeviceToDevice);
-
-	gT_FanVolumeTexture.addressMode[0] = cudaAddressModeBorder;
-	gT_FanVolumeTexture.addressMode[1] = cudaAddressModeBorder;
-	gT_FanVolumeTexture.filterMode = cudaFilterModeLinear;
-	gT_FanVolumeTexture.normalized = false;
-
-	// TODO: For very small sizes (roughly <=512x128) with few angles (<=180)
-	// not using an array is more efficient.
-	//cudaBindTexture2D(0, gT_FanVolumeTexture, (const void*)data, channelDesc, width, height, sizeof(float)*pitch);
-	cudaBindTextureToArray(gT_FanVolumeTexture, dataArray, channelDesc);
-
-	// TODO: error value?
-
-	return true;
-}
-
 // projection for angles that are roughly horizontal
 // (detector roughly vertical)
-__global__ void FanFPhorizontal(float* D_projData, unsigned int projPitch, unsigned int startSlice, unsigned int startAngle, unsigned int endAngle, const SDimensions dims, float outputScale)
+__global__ void FanFPhorizontal(float* D_projData, unsigned int projPitch, cudaTextureObject_t tex, unsigned int startSlice, unsigned int startAngle, unsigned int endAngle, const SDimensions dims, float outputScale)
 {
 	float* projData = (float*)D_projData;
 	const int relDet = threadIdx.x;
@@ -134,7 +106,7 @@ __global__ void FanFPhorizontal(float* D_projData, unsigned int projPitch, unsig
 		float fV = 0.0f;
 		for (int slice = startSlice; slice < endSlice; ++slice)
 		{
-			fV += tex2D(gT_FanVolumeTexture, fX, fY);
+			fV += tex2D<float>(tex, fX, fY);
 			fY -= alpha;
 			fX += 1.0f;
 		}
@@ -149,7 +121,7 @@ __global__ void FanFPhorizontal(float* D_projData, unsigned int projPitch, unsig
 
 // projection for angles that are roughly vertical
 // (detector roughly horizontal)
-__global__ void FanFPvertical(float* D_projData, unsigned int projPitch, unsigned int startSlice, unsigned int startAngle, unsigned int endAngle, const SDimensions dims, float outputScale)
+__global__ void FanFPvertical(float* D_projData, unsigned int projPitch, cudaTextureObject_t tex, unsigned int startSlice, unsigned int startAngle, unsigned int endAngle, const SDimensions dims, float outputScale)
 {
 	const int relDet = threadIdx.x;
 	const int relAngle = threadIdx.y;
@@ -207,7 +179,7 @@ __global__ void FanFPvertical(float* D_projData, unsigned int projPitch, unsigne
 
 		for (int slice = startSlice; slice < endSlice; ++slice)
 		{
-			fV += tex2D(gT_FanVolumeTexture, fX, fY);
+			fV += tex2D<float>(tex, fX, fY);
 			fX -= alpha;
 			fY += 1.0f;
 		}
@@ -227,7 +199,10 @@ bool FanFP_internal(float* D_volumeData, unsigned int volumePitch,
 	assert(dims.iProjAngles <= g_MaxAngles);
 
 	cudaArray* D_dataArray;
-	bindVolumeDataTexture(D_volumeData, D_dataArray, volumePitch, dims.iVolWidth, dims.iVolHeight);
+	cudaTextureObject_t D_texObj;
+
+	if (!createArrayAndTextureObject2D(D_volumeData, D_dataArray, D_texObj, volumePitch, dims.iVolWidth, dims.iVolHeight))
+		return false;
 
 	// transfer angles to constant memory
 	float* tmp = new float[dims.iProjAngles];
@@ -260,13 +235,13 @@ bool FanFP_internal(float* D_volumeData, unsigned int volumePitch,
 	cudaStreamCreate(&stream1);
 	streams.push_back(stream1);
 	for (unsigned int i = 0; i < dims.iVolWidth; i += g_blockSlices)
-		FanFPhorizontal<<<dimGrid, dimBlock, 0, stream1>>>(D_projData, projPitch, i, blockStart, blockEnd, dims, outputScale);
+		FanFPhorizontal<<<dimGrid, dimBlock, 0, stream1>>>(D_projData, projPitch, D_texObj, i, blockStart, blockEnd, dims, outputScale);
 
 	cudaStream_t stream2;
 	cudaStreamCreate(&stream2);
 	streams.push_back(stream2);
 	for (unsigned int i = 0; i < dims.iVolHeight; i += g_blockSlices)
-		FanFPvertical<<<dimGrid, dimBlock, 0, stream2>>>(D_projData, projPitch, i, blockStart, blockEnd, dims, outputScale);
+		FanFPvertical<<<dimGrid, dimBlock, 0, stream2>>>(D_projData, projPitch, D_texObj, i, blockStart, blockEnd, dims, outputScale);
 
 	bool ok = true;
 
@@ -278,6 +253,8 @@ bool FanFP_internal(float* D_volumeData, unsigned int volumePitch,
 
 	cudaFreeArray(D_dataArray);
 
+	cudaDestroyTextureObject(D_texObj);
+
 	return ok;
 }
 
diff --git a/cuda/2d/par_bp.cu b/cuda/2d/par_bp.cu
index d7c3ab0..4066912 100644
--- a/cuda/2d/par_bp.cu
+++ b/cuda/2d/par_bp.cu
@@ -51,24 +51,8 @@ __constant__ float gC_angle_scaled_cos[g_MaxAngles];
 __constant__ float gC_angle_offset[g_MaxAngles];
 __constant__ float gC_angle_scale[g_MaxAngles];
 
-static bool bindProjDataTexture(float* data, unsigned int pitch, unsigned int width, unsigned int height, cudaTextureAddressMode mode = cudaAddressModeBorder)
-{
-	cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<float>();
-
-	gT_projTexture.addressMode[0] = mode;
-	gT_projTexture.addressMode[1] = mode;
-	gT_projTexture.filterMode = cudaFilterModeLinear;
-	gT_projTexture.normalized = false;
-
-	cudaBindTexture2D(0, gT_projTexture, (const void*)data, channelDesc, width, height, sizeof(float)*pitch);
-
-	// TODO: error value?
-
-	return true;
-}
-
 // TODO: Templated version with/without scale? (Or only the global outputscale)
-__global__ void devBP(float* D_volData, unsigned int volPitch, unsigned int startAngle, const SDimensions dims, float fOutputScale)
+__global__ void devBP(float* D_volData, unsigned int volPitch, cudaTextureObject_t tex, unsigned int startAngle, const SDimensions dims, float fOutputScale)
 {
 	const int relX = threadIdx.x;
 	const int relY = threadIdx.y;
@@ -98,7 +82,7 @@ __global__ void devBP(float* D_volData, unsigned int volPitch, unsigned int star
 		const float scale = gC_angle_scale[angle];
 
 		const float fT = fX * scaled_cos_theta - fY * scaled_sin_theta + TOffset;
-		fVal += tex2D(gT_projTexture, fT, fA) * scale;
+		fVal += tex2D<float>(tex, fT, fA) * scale;
 		fA += 1.0f;
 	}
 
@@ -106,7 +90,7 @@ __global__ void devBP(float* D_volData, unsigned int volPitch, unsigned int star
 }
 
 // supersampling version
-__global__ void devBP_SS(float* D_volData, unsigned int volPitch, unsigned int startAngle, const SDimensions dims, float fOutputScale)
+__global__ void devBP_SS(float* D_volData, unsigned int volPitch, cudaTextureObject_t tex, unsigned int startAngle, const SDimensions dims, float fOutputScale)
 {
 	const int relX = threadIdx.x;
 	const int relY = threadIdx.y;
@@ -145,7 +129,7 @@ __global__ void devBP_SS(float* D_volData, unsigned int volPitch, unsigned int s
 			float fTy = fT;
 			fT += fSubStep * cos_theta;
 			for (int iSubY = 0; iSubY < dims.iRaysPerPixelDim; ++iSubY) {
-				fVal += tex2D(gT_projTexture, fTy, fA) * scale;
+				fVal += tex2D<float>(tex, fTy, fA) * scale;
 				fTy -= fSubStep * sin_theta;
 			}
 		}
@@ -155,7 +139,7 @@ __global__ void devBP_SS(float* D_volData, unsigned int volPitch, unsigned int s
 	volData[Y*volPitch+X] += fVal * fOutputScale;
 }
 
-__global__ void devBP_SART(float* D_volData, unsigned int volPitch, float offset, float angle_sin, float angle_cos, const SDimensions dims, float fOutputScale)
+__global__ void devBP_SART(float* D_volData, unsigned int volPitch, cudaTextureObject_t tex, float offset, float angle_sin, float angle_cos, const SDimensions dims, float fOutputScale)
 {
 	const int relX = threadIdx.x;
 	const int relY = threadIdx.y;
@@ -170,7 +154,7 @@ __global__ void devBP_SART(float* D_volData, unsigned int volPitch, float offset
 	const float fY = ( Y - 0.5f*dims.iVolHeight + 0.5f );
 
 	const float fT = fX * angle_cos - fY * angle_sin + offset;
-	const float fVal = tex2D(gT_projTexture, fT, 0.5f);
+	const float fVal = tex2D<float>(tex, fT, 0.5f);
 
 	// NB: The 'scale' constant in devBP is cancelled out by the SART weighting
 
@@ -185,13 +169,15 @@ bool BP_internal(float* D_volumeData, unsigned int volumePitch,
 {
 	assert(dims.iProjAngles <= g_MaxAngles);
 
+	cudaTextureObject_t D_texObj;
+	if (!createTextureObjectPitch2D(D_projData, D_texObj, projPitch, dims.iProjDets, dims.iProjAngles))
+		return false;
+
 	float* angle_scaled_sin = new float[dims.iProjAngles];
 	float* angle_scaled_cos = new float[dims.iProjAngles];
 	float* angle_offset = new float[dims.iProjAngles];
 	float* angle_scale = new float[dims.iProjAngles];
 
-	bindProjDataTexture(D_projData, projPitch, dims.iProjDets, dims.iProjAngles);
-
 	for (unsigned int i = 0; i < dims.iProjAngles; ++i) {
 		double d = angles[i].fDetUX * angles[i].fRayY - angles[i].fDetUY * angles[i].fRayX;
 		angle_scaled_cos[i] = angles[i].fRayY / d;
@@ -227,15 +213,17 @@ bool BP_internal(float* D_volumeData, unsigned int volumePitch,
 	for (unsigned int i = 0; i < dims.iProjAngles; i += g_anglesPerBlock) {
 
 		if (dims.iRaysPerPixelDim > 1)
-			devBP_SS<<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, i, dims, fOutputScale);
+			devBP_SS<<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, D_texObj, i, dims, fOutputScale);
 		else
-			devBP<<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, i, dims, fOutputScale);
+			devBP<<<dimGrid, dimBlock, 0, stream>>>(D_volumeData, volumePitch, D_texObj, i, dims, fOutputScale);
 	}
 
 	bool ok = checkCuda(cudaStreamSynchronize(stream), "par_bp");
 
 	cudaStreamDestroy(stream);
 
+	cudaDestroyTextureObject(D_texObj);
+
 	return ok;
 }
 
@@ -269,7 +257,9 @@ bool BP_SART(float* D_volumeData, unsigned int volumePitch,
 	// Only one angle.
 	// We need to Clamp to the border pixels instead of to zero, because
 	// SART weights with ray length.
-	bindProjDataTexture(D_projData, projPitch, dims.iProjDets, 1, cudaAddressModeClamp);
+	cudaTextureObject_t D_texObj;
+	if (!createTextureObjectPitch2D(D_projData, D_texObj, projPitch, dims.iProjDets, 1, cudaAddressModeClamp))
+		return false;
 
 	double d = angles[angle].fDetUX * angles[angle].fRayY - angles[angle].fDetUY * angles[angle].fRayX;
 	float angle_scaled_cos = angles[angle].fRayY / d;
@@ -282,9 +272,13 @@ bool BP_SART(float* D_volumeData, unsigned int volumePitch,
 	dim3 dimGrid((dims.iVolWidth+g_blockSlices-1)/g_blockSlices,
 	             (dims.iVolHeight+g_blockSliceSize-1)/g_blockSliceSize);
 
-	devBP_SART<<<dimGrid, dimBlock>>>(D_volumeData, volumePitch, angle_offset, angle_scaled_sin, angle_scaled_cos, dims, fOutputScale);
+	devBP_SART<<<dimGrid, dimBlock>>>(D_volumeData, volumePitch, D_texObj, angle_offset, angle_scaled_sin, angle_scaled_cos, dims, fOutputScale);
 
-	return checkCuda(cudaThreadSynchronize(), "BP_SART");
+	bool ok = checkCuda(cudaThreadSynchronize(), "BP_SART");
+
+	cudaDestroyTextureObject(D_texObj);
+
+	return ok;
 }
 
 
diff --git a/cuda/2d/par_fp.cu b/cuda/2d/par_fp.cu
index e947428..3d9e842 100644
--- a/cuda/2d/par_fp.cu
+++ b/cuda/2d/par_fp.cu
@@ -34,11 +34,6 @@ along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
 #include <list>
 #include <cmath>
 
-typedef texture<float, 2, cudaReadModeElementType> texture2D;
-
-static texture2D gT_volumeTexture;
-
-
 namespace astraCUDA {
 
 static const unsigned g_MaxAngles = 2560;
@@ -52,37 +47,9 @@ static const unsigned int g_anglesPerBlock = 16;
 static const unsigned int g_detBlockSize = 32;
 static const unsigned int g_blockSlices = 64;
 
-// fixed point scaling factor
-#define fPREC_FACTOR 16.0f
-#define iPREC_FACTOR 16
-
-
-static bool bindVolumeDataTexture(float* data, cudaArray*& dataArray, unsigned int pitch, unsigned int width, unsigned int height)
-{
-	cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<float>();
-	dataArray = 0;
-	cudaMallocArray(&dataArray, &channelDesc, width, height);
-	cudaMemcpy2DToArray(dataArray, 0, 0, data, pitch*sizeof(float), width*sizeof(float), height, cudaMemcpyDeviceToDevice);
-
-	gT_volumeTexture.addressMode[0] = cudaAddressModeBorder;
-	gT_volumeTexture.addressMode[1] = cudaAddressModeBorder;
-	gT_volumeTexture.filterMode = cudaFilterModeLinear;
-	gT_volumeTexture.normalized = false;
-
-	// TODO: For very small sizes (roughly <=512x128) with few angles (<=180)
-	// not using an array is more efficient.
-//	cudaBindTexture2D(0, gT_volumeTexture, (const void*)data, channelDesc, width, height, sizeof(float)*pitch);
-	cudaBindTextureToArray(gT_volumeTexture, dataArray, channelDesc);
-
-	// TODO: error value?
-
-	return true;
-}
-
-
 // projection for angles that are roughly horizontal
 // (detector roughly vertical)
-__global__ void FPhorizontal_simple(float* D_projData, unsigned int projPitch, unsigned int startSlice, unsigned int startAngle, unsigned int endAngle, const SDimensions dims, float outputScale)
+__global__ void FPhorizontal_simple(float* D_projData, unsigned int projPitch, cudaTextureObject_t tex, unsigned int startSlice, unsigned int startAngle, unsigned int endAngle, const SDimensions dims, float outputScale)
 {
 	const int relDet = threadIdx.x;
 	const int relAngle = threadIdx.y;
@@ -134,7 +101,7 @@ __global__ void FPhorizontal_simple(float* D_projData, unsigned int projPitch, u
 		for (int slice = startSlice; slice < endSlice; ++slice)
 		{
 			for (int iSubT = 0; iSubT < dims.iRaysPerDet; ++iSubT) {
-				fVal += tex2D(gT_volumeTexture, fS, fP);
+				fVal += tex2D<float>(tex, fS, fP);
 				fP += fSubDetStep;
 			}
 			fP += fSliceStep;
@@ -145,7 +112,7 @@ __global__ void FPhorizontal_simple(float* D_projData, unsigned int projPitch, u
 
 		for (int slice = startSlice; slice < endSlice; ++slice)
 		{
-			fVal += tex2D(gT_volumeTexture, fS, fP);
+			fVal += tex2D<float>(tex, fS, fP);
 			fP += fSliceStep;
 			fS += 1.0f;
 		}
@@ -159,7 +126,7 @@ __global__ void FPhorizontal_simple(float* D_projData, unsigned int projPitch, u
 
 // projection for angles that are roughly vertical
 // (detector roughly horizontal)
-__global__ void FPvertical_simple(float* D_projData, unsigned int projPitch, unsigned int startSlice, unsigned int startAngle, unsigned int endAngle, const SDimensions dims, float outputScale)
+__global__ void FPvertical_simple(float* D_projData, unsigned int projPitch, cudaTextureObject_t tex, unsigned int startSlice, unsigned int startAngle, unsigned int endAngle, const SDimensions dims, float outputScale)
 {
 	const int relDet = threadIdx.x;
 	const int relAngle = threadIdx.y;
@@ -210,7 +177,7 @@ __global__ void FPvertical_simple(float* D_projData, unsigned int projPitch, uns
 		for (int slice = startSlice; slice < endSlice; ++slice)
 		{
 			for (int iSubT = 0; iSubT < dims.iRaysPerDet; ++iSubT) {
-				fVal += tex2D(gT_volumeTexture, fP, fS);
+				fVal += tex2D<float>(tex, fP, fS);
 				fP += fSubDetStep;
 			}
 			fP += fSliceStep;
@@ -221,7 +188,7 @@ __global__ void FPvertical_simple(float* D_projData, unsigned int projPitch, uns
 
 		for (int slice = startSlice; slice < endSlice; ++slice)
 		{
-			fVal += tex2D(gT_volumeTexture, fP, fS);
+			fVal += tex2D<float>(tex, fP, fS);
 			fP += fSliceStep;
 			fS += 1.0f;
 		}
@@ -269,7 +236,10 @@ bool FP_simple_internal(float* D_volumeData, unsigned int volumePitch,
 	assert(angles);
 
 	cudaArray* D_dataArray;
-	bindVolumeDataTexture(D_volumeData, D_dataArray, volumePitch, dims.iVolWidth, dims.iVolHeight);
+	cudaTextureObject_t D_texObj;
+
+	if (!createArrayAndTextureObject2D(D_volumeData, D_dataArray, D_texObj, volumePitch, dims.iVolWidth, dims.iVolHeight))
+		return false;
 
 
 	convertAndUploadAngles(angles, dims.iProjAngles, dims.iProjDets);
@@ -313,10 +283,10 @@ bool FP_simple_internal(float* D_volumeData, unsigned int volumePitch,
 				//printf("angle block: %d to %d, %d\n", blockStart, blockEnd, blockVertical);
 				if (!blockVertical)
 					for (unsigned int i = 0; i < dims.iVolWidth; i += g_blockSlices)
-						FPhorizontal_simple<<<dimGrid, dimBlock, 0, stream>>>(D_projData, projPitch, i, blockStart, blockEnd, dims, outputScale);
+						FPhorizontal_simple<<<dimGrid, dimBlock, 0, stream>>>(D_projData, projPitch, D_texObj, i, blockStart, blockEnd, dims, outputScale);
 				else
 					for (unsigned int i = 0; i < dims.iVolHeight; i += g_blockSlices)
-						FPvertical_simple<<<dimGrid, dimBlock, 0, stream>>>(D_projData, projPitch, i, blockStart, blockEnd, dims, outputScale);
+						FPvertical_simple<<<dimGrid, dimBlock, 0, stream>>>(D_projData, projPitch, D_texObj, i, blockStart, blockEnd, dims, outputScale);
 			}
 			blockVertical = vertical;
 			blockStart = a;
@@ -332,6 +302,8 @@ bool FP_simple_internal(float* D_volumeData, unsigned int volumePitch,
 
 	cudaFreeArray(D_dataArray);
 
+	cudaDestroyTextureObject(D_texObj);
+
 	return ok;
 }
 
diff --git a/cuda/2d/util.cu b/cuda/2d/util.cu
index ac360f0..4a58880 100644
--- a/cuda/2d/util.cu
+++ b/cuda/2d/util.cu
@@ -126,6 +126,75 @@ void duplicateProjectionData(float* D_dst, float* D_src, unsigned int pitch, con
 	cudaMemcpy2D(D_dst, sizeof(float)*pitch, D_src, sizeof(float)*pitch, sizeof(float)*dims.iProjDets, dims.iProjAngles, cudaMemcpyDeviceToDevice);
 }
 
+bool createArrayAndTextureObject2D(float* data, cudaArray*& dataArray, cudaTextureObject_t& texObj, unsigned int pitch, unsigned int width, unsigned int height)
+{
+	// TODO: For very small sizes (roughly <=512x128) with few angles (<=180)
+	// not using an array is more efficient.
+
+	cudaChannelFormatDesc channelDesc =
+	    cudaCreateChannelDesc(32, 0, 0, 0, cudaChannelFormatKindFloat);
+
+	dataArray = 0;
+	if (!checkCuda(cudaMallocArray(&dataArray, &channelDesc, width, height), "createTextureObject2D malloc"))
+		return false;
+	if (!checkCuda(cudaMemcpy2DToArray(dataArray, 0, 0, data, pitch*sizeof(float), width*sizeof(float), height, cudaMemcpyDeviceToDevice), "createTextureObject2D memcpy")) {
+		cudaFreeArray(dataArray);
+		return false;
+	}
+
+	cudaResourceDesc resDesc;
+	memset(&resDesc, 0, sizeof(resDesc));
+	resDesc.resType = cudaResourceTypeArray;
+	resDesc.res.array.array = dataArray;
+
+	cudaTextureDesc texDesc;
+	memset(&texDesc, 0, sizeof(texDesc));
+	texDesc.addressMode[0] = cudaAddressModeBorder;
+	texDesc.addressMode[1] = cudaAddressModeBorder;
+	texDesc.filterMode = cudaFilterModeLinear;
+	texDesc.readMode = cudaReadModeElementType;
+	texDesc.normalizedCoords = 0;
+
+	texObj = 0;
+
+	if (!checkCuda(cudaCreateTextureObject(&texObj, &resDesc, &texDesc, NULL), "createTextureObject2D")) {
+		cudaFreeArray(dataArray);
+		return false;
+	}
+
+	return true;
+}
+
+bool createTextureObjectPitch2D(float* data, cudaTextureObject_t& texObj, unsigned int pitch, unsigned int width, unsigned int height, cudaTextureAddressMode mode)
+{
+	cudaChannelFormatDesc channelDesc =
+	    cudaCreateChannelDesc(32, 0, 0, 0, cudaChannelFormatKindFloat);
+
+	cudaResourceDesc resDesc;
+	memset(&resDesc, 0, sizeof(resDesc));
+	resDesc.resType = cudaResourceTypePitch2D;
+	resDesc.res.pitch2D.devPtr = (void*)data;
+	resDesc.res.pitch2D.desc = channelDesc;
+	resDesc.res.pitch2D.width = width;
+	resDesc.res.pitch2D.height = height;
+	resDesc.res.pitch2D.pitchInBytes = sizeof(float)*pitch;
+
+	cudaTextureDesc texDesc;
+	memset(&texDesc, 0, sizeof(texDesc));
+	texDesc.addressMode[0] = mode;
+	texDesc.addressMode[1] = mode;
+	texDesc.filterMode = cudaFilterModeLinear;
+	texDesc.readMode = cudaReadModeElementType;
+	texDesc.normalizedCoords = 0;
+
+	texObj = 0;
+
+	return checkCuda(cudaCreateTextureObject(&texObj, &resDesc, &texDesc, NULL), "createTextureObjectPitch2D");
+}
+
+
+
+
 template <unsigned int blockSize>
 __global__ void reduce1D(float *g_idata, float *g_odata, unsigned int n)
 {
diff --git a/cuda/3d/cone_bp.cu b/cuda/3d/cone_bp.cu
index e265304..41d781c 100644
--- a/cuda/3d/cone_bp.cu
+++ b/cuda/3d/cone_bp.cu
@@ -35,10 +35,6 @@ along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
 
 #include <cuda.h>
 
-typedef texture<float, 3, cudaReadModeElementType> texture3D;
-
-static texture3D gT_coneProjTexture;
-
 namespace astraCUDA3d {
 
 static const unsigned int g_volBlockZ = 6;
@@ -57,28 +53,12 @@ struct DevConeParams {
 
 __constant__ DevConeParams gC_C[g_MaxAngles];
 
-bool bindProjDataTexture(const cudaArray* array)
-{
-	cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<float>();
-
-	gT_coneProjTexture.addressMode[0] = cudaAddressModeBorder;
-	gT_coneProjTexture.addressMode[1] = cudaAddressModeBorder;
-	gT_coneProjTexture.addressMode[2] = cudaAddressModeBorder;
-	gT_coneProjTexture.filterMode = cudaFilterModeLinear;
-	gT_coneProjTexture.normalized = false;
-
-	cudaBindTextureToArray(gT_coneProjTexture, array, channelDesc);
-
-	// TODO: error value?
-
-	return true;
-}
-
-
 //__launch_bounds__(32*16, 4)
 template<bool FDKWEIGHT, unsigned int ZSIZE>
-__global__ void dev_cone_BP(void* D_volData, unsigned int volPitch, int startAngle,
-                            int angleOffset, const astraCUDA3d::SDimensions3D dims,
+__global__ void dev_cone_BP(void* D_volData, unsigned int volPitch,
+                            cudaTextureObject_t tex,
+                            int startAngle, int angleOffset,
+                            const astraCUDA3d::SDimensions3D dims,
                             float fOutputScale)
 {
 	float* volData = (float*)D_volData;
@@ -133,7 +113,7 @@ __global__ void dev_cone_BP(void* D_volData, unsigned int volPitch, int startAng
 				fr = __fdividef(1.0f, fDen);
 				fU = fUNum * fr;
 				fV = fVNum * fr;
-				float fVal = tex3D(gT_coneProjTexture, fU, fAngle, fV);
+				float fVal = tex3D<float>(tex, fU, fAngle, fV);
 				Z[idx] += fr*fr*fVal;
 
 				fUNum += fCu.z;
@@ -154,7 +134,7 @@ __global__ void dev_cone_BP(void* D_volData, unsigned int volPitch, int startAng
 
 
 // supersampling version
-__global__ void dev_cone_BP_SS(void* D_volData, unsigned int volPitch, int startAngle, int angleOffset, const SDimensions3D dims, int iRaysPerVoxelDim, float fOutputScale)
+__global__ void dev_cone_BP_SS(void* D_volData, unsigned int volPitch, cudaTextureObject_t tex, int startAngle, int angleOffset, const SDimensions3D dims, int iRaysPerVoxelDim, float fOutputScale)
 {
 	float* volData = (float*)D_volData;
 
@@ -220,7 +200,7 @@ __global__ void dev_cone_BP_SS(void* D_volData, unsigned int volPitch, int start
 				const float fU = fUNum * fr;
 				const float fV = fVNum * fr;
 
-				fVal += tex3D(gT_coneProjTexture, fU, fAngle, fV) * fr * fr;
+				fVal += tex3D<float>(tex, fU, fAngle, fV) * fr * fr;
 
 				fZs += fSubStep;
 			}
@@ -313,7 +293,9 @@ bool ConeBP_Array(cudaPitchedPtr D_volumeData,
                   const SDimensions3D& dims, const SConeProjection* angles,
                   const SProjectorParams3D& params)
 {
-	bindProjDataTexture(D_projArray);
+	cudaTextureObject_t D_texObj;
+	if (!createTextureObject3D(D_projArray, D_texObj))
+		return false;
 
 	float fOutputScale;
 	if (params.bFDKWeighting) {
@@ -323,14 +305,16 @@ bool ConeBP_Array(cudaPitchedPtr D_volumeData,
 		fOutputScale = params.fOutputScale * (params.fVolScaleX * params.fVolScaleY * params.fVolScaleZ);
 	}
 
+	bool ok = true;
+
 	for (unsigned int th = 0; th < dims.iProjAngles; th += g_MaxAngles) {
 		unsigned int angleCount = g_MaxAngles;
 		if (th + angleCount > dims.iProjAngles)
 			angleCount = dims.iProjAngles - th;
 
-		bool ok = transferConstants(angles, angleCount, params);
+		ok = transferConstants(angles, angleCount, params);
 		if (!ok)
-			return false;
+			break;
 
 		dim3 dimBlock(g_volBlockX, g_volBlockY);
 
@@ -343,31 +327,33 @@ bool ConeBP_Array(cudaPitchedPtr D_volumeData,
 		// printf("Calling BP: %d, %dx%d, %dx%d to %p\n", i, dimBlock.x, dimBlock.y, dimGrid.x, dimGrid.y, (void*)D_volumeData.ptr); 
 			if (params.bFDKWeighting) {
 				if (dims.iVolZ == 1) {
-					dev_cone_BP<true, 1><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), i, th, dims, fOutputScale);
+					dev_cone_BP<true, 1><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), D_texObj, i, th, dims, fOutputScale);
 				} else {
-					dev_cone_BP<true, g_volBlockZ><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), i, th, dims, fOutputScale);
+					dev_cone_BP<true, g_volBlockZ><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), D_texObj, i, th, dims, fOutputScale);
 				}
 			} else if (params.iRaysPerVoxelDim == 1) {
 				if (dims.iVolZ == 1) {
-					dev_cone_BP<false, 1><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), i, th, dims, fOutputScale);
+					dev_cone_BP<false, 1><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), D_texObj, i, th, dims, fOutputScale);
 				} else {
-					dev_cone_BP<false, g_volBlockZ><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), i, th, dims, fOutputScale);
+					dev_cone_BP<false, g_volBlockZ><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), D_texObj, i, th, dims, fOutputScale);
 				}
 			} else
-				dev_cone_BP_SS<<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), i, th, dims, params.iRaysPerVoxelDim, fOutputScale);
+				dev_cone_BP_SS<<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), D_texObj, i, th, dims, params.iRaysPerVoxelDim, fOutputScale);
 		}
 
 		// TODO: Consider not synchronizing here, if possible.
-		if (!checkCuda(cudaThreadSynchronize(), "cone_bp"))
-			return false;
+		ok = checkCuda(cudaThreadSynchronize(), "cone_bp");
+		if (!ok)
+			break;
 
 		angles = angles + angleCount;
 		// printf("%f\n", toc(t));
 
 	}
 
+	cudaDestroyTextureObject(D_texObj);
 
-	return true;
+	return ok;
 }
 
 bool ConeBP(cudaPitchedPtr D_volumeData,
diff --git a/cuda/3d/cone_fp.cu b/cuda/3d/cone_fp.cu
index 2c3d1f6..2ef58ee 100644
--- a/cuda/3d/cone_fp.cu
+++ b/cuda/3d/cone_fp.cu
@@ -35,10 +35,6 @@ along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
 
 #include <cuda.h>
 
-typedef texture<float, 3, cudaReadModeElementType> texture3D;
-
-static texture3D gT_coneVolumeTexture;
-
 namespace astraCUDA3d {
 
 static const unsigned int g_anglesPerBlock = 4;
@@ -63,24 +59,6 @@ __constant__ float gC_DetVY[g_MaxAngles];
 __constant__ float gC_DetVZ[g_MaxAngles];
 
 
-bool bindVolumeDataTexture(const cudaArray* array)
-{
-	cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<float>();
-
-	gT_coneVolumeTexture.addressMode[0] = cudaAddressModeBorder;
-	gT_coneVolumeTexture.addressMode[1] = cudaAddressModeBorder;
-	gT_coneVolumeTexture.addressMode[2] = cudaAddressModeBorder;
-	gT_coneVolumeTexture.filterMode = cudaFilterModeLinear;
-	gT_coneVolumeTexture.normalized = false;
-
-	cudaBindTextureToArray(gT_coneVolumeTexture, array, channelDesc);
-
-	// TODO: error value?
-
-	return true;
-}
-
-
 // x=0, y=1, z=2
 struct DIR_X {
 	__device__ float nSlices(const SDimensions3D& dims) const { return dims.iVolX; }
@@ -89,7 +67,7 @@ struct DIR_X {
 	__device__ float c0(float x, float y, float z) const { return x; }
 	__device__ float c1(float x, float y, float z) const { return y; }
 	__device__ float c2(float x, float y, float z) const { return z; }
-	__device__ float tex(float f0, float f1, float f2) const { return tex3D(gT_coneVolumeTexture, f0, f1, f2); }
+	__device__ float tex(cudaTextureObject_t tex, float f0, float f1, float f2) const { return tex3D<float>(tex, f0, f1, f2); }
 	__device__ float x(float f0, float f1, float f2) const { return f0; }
 	__device__ float y(float f0, float f1, float f2) const { return f1; }
 	__device__ float z(float f0, float f1, float f2) const { return f2; }
@@ -103,7 +81,7 @@ struct DIR_Y {
 	__device__ float c0(float x, float y, float z) const { return y; }
 	__device__ float c1(float x, float y, float z) const { return x; }
 	__device__ float c2(float x, float y, float z) const { return z; }
-	__device__ float tex(float f0, float f1, float f2) const { return tex3D(gT_coneVolumeTexture, f1, f0, f2); }
+	__device__ float tex(cudaTextureObject_t tex, float f0, float f1, float f2) const { return tex3D<float>(tex, f1, f0, f2); }
 	__device__ float x(float f0, float f1, float f2) const { return f1; }
 	__device__ float y(float f0, float f1, float f2) const { return f0; }
 	__device__ float z(float f0, float f1, float f2) const { return f2; }
@@ -117,7 +95,7 @@ struct DIR_Z {
 	__device__ float c0(float x, float y, float z) const { return z; }
 	__device__ float c1(float x, float y, float z) const { return x; }
 	__device__ float c2(float x, float y, float z) const { return y; }
-	__device__ float tex(float f0, float f1, float f2) const { return tex3D(gT_coneVolumeTexture, f1, f2, f0); }
+	__device__ float tex(cudaTextureObject_t tex, float f0, float f1, float f2) const { return tex3D<float>(tex, f1, f2, f0); }
 	__device__ float x(float f0, float f1, float f2) const { return f1; }
 	__device__ float y(float f0, float f1, float f2) const { return f2; }
 	__device__ float z(float f0, float f1, float f2) const { return f0; }
@@ -144,6 +122,7 @@ struct SCALE_NONCUBE {
 
 template<class COORD, class SCALE>
 __global__ void cone_FP_t(float* D_projData, unsigned int projPitch,
+                          cudaTextureObject_t tex,
                           unsigned int startSlice,
                           unsigned int startAngle, unsigned int endAngle,
                           const SDimensions3D dims,
@@ -208,7 +187,7 @@ __global__ void cone_FP_t(float* D_projData, unsigned int projPitch,
 
 		for (int s = startSlice; s < endSlice; ++s)
 		{
-			fVal += c.tex(f0, f1, f2);
+			fVal += c.tex(tex, f0, f1, f2);
 			f0 += 1.0f;
 			f1 += a1;
 			f2 += a2;
@@ -222,6 +201,7 @@ __global__ void cone_FP_t(float* D_projData, unsigned int projPitch,
 
 template<class COORD>
 __global__ void cone_FP_SS_t(float* D_projData, unsigned int projPitch,
+                             cudaTextureObject_t tex,
                              unsigned int startSlice,
                              unsigned int startAngle, unsigned int endAngle,
                              const SDimensions3D dims, int iRaysPerDetDim,
@@ -295,7 +275,7 @@ __global__ void cone_FP_SS_t(float* D_projData, unsigned int projPitch,
 
 		for (int s = startSlice; s < endSlice; ++s)
 		{
-			fVal += c.tex(f0, f1, f2);
+			fVal += c.tex(tex, f0, f1, f2);
 			f0 += 1.0f;
 			f1 += a1;
 			f2 += a2;
@@ -313,7 +293,9 @@ __global__ void cone_FP_SS_t(float* D_projData, unsigned int projPitch,
 
 
 bool ConeFP_Array_internal(cudaPitchedPtr D_projData,
-                  const SDimensions3D& dims, unsigned int angleCount, const SConeProjection* angles,
+                  cudaTextureObject_t D_texObj,
+                  const SDimensions3D& dims,
+                  unsigned int angleCount, const SConeProjection* angles,
                   const SProjectorParams3D& params)
 {
 	// transfer angles to constant memory
@@ -419,29 +401,29 @@ bool ConeFP_Array_internal(cudaPitchedPtr D_projData,
 					for (unsigned int i = 0; i < dims.iVolX; i += g_blockSlices)
 						if (params.iRaysPerDetDim == 1)
 							if (cube)
-								cone_FP_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, scube);
+								cone_FP_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, scube);
 							else
-								cone_FP_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, snoncubeX);
+								cone_FP_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, snoncubeX);
 						else
-							cone_FP_SS_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeX);
+							cone_FP_SS_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeX);
 				} else if (blockDirection == 1) {
 					for (unsigned int i = 0; i < dims.iVolY; i += g_blockSlices)
 						if (params.iRaysPerDetDim == 1)
 							if (cube)
-								cone_FP_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, scube);
+								cone_FP_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, scube);
 							else
-								cone_FP_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, snoncubeY);
+								cone_FP_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, snoncubeY);
 						else
-							cone_FP_SS_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeY);
+							cone_FP_SS_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeY);
 				} else if (blockDirection == 2) {
 					for (unsigned int i = 0; i < dims.iVolZ; i += g_blockSlices)
 						if (params.iRaysPerDetDim == 1)
 							if (cube)
-								cone_FP_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, scube);
+								cone_FP_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, scube);
 							else
-								cone_FP_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, snoncubeZ);
+								cone_FP_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, snoncubeZ);
 						else
-							cone_FP_SS_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeZ);
+							cone_FP_SS_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeZ);
 				}
 
 			}
@@ -470,10 +452,14 @@ bool ConeFP(cudaPitchedPtr D_volumeData,
             const SProjectorParams3D& params)
 {
 	// transfer volume to array
-
 	cudaArray* cuArray = allocateVolumeArray(dims);
 	transferVolumeToArray(D_volumeData, cuArray, dims);
-	bindVolumeDataTexture(cuArray);
+
+	cudaTextureObject_t D_texObj;
+	if (!createTextureObject3D(cuArray, D_texObj)) {
+		cudaFreeArray(cuArray);
+		return false;
+	}
 
 	bool ret;
 
@@ -485,7 +471,7 @@ bool ConeFP(cudaPitchedPtr D_volumeData,
 		cudaPitchedPtr D_subprojData = D_projData;
 		D_subprojData.ptr = (char*)D_projData.ptr + iAngle * D_projData.pitch;
 
-		ret = ConeFP_Array_internal(D_subprojData,
+		ret = ConeFP_Array_internal(D_subprojData, D_texObj,
 		                            dims, iEndAngle - iAngle, angles + iAngle,
 		                            params);
 		if (!ret)
diff --git a/cuda/3d/par3d_bp.cu b/cuda/3d/par3d_bp.cu
index 1dc75ce..27d95fe 100644
--- a/cuda/3d/par3d_bp.cu
+++ b/cuda/3d/par3d_bp.cu
@@ -35,10 +35,6 @@ along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
 
 #include <cuda.h>
 
-typedef texture<float, 3, cudaReadModeElementType> texture3D;
-
-static texture3D gT_par3DProjTexture;
-
 namespace astraCUDA3d {
 
 static const unsigned int g_volBlockZ = 6;
@@ -58,26 +54,8 @@ __constant__ DevPar3DParams gC_C[g_MaxAngles];
 __constant__ float gC_scale[g_MaxAngles];
 
 
-static bool bindProjDataTexture(const cudaArray* array)
-{
-	cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<float>();
-
-	gT_par3DProjTexture.addressMode[0] = cudaAddressModeBorder;
-	gT_par3DProjTexture.addressMode[1] = cudaAddressModeBorder;
-	gT_par3DProjTexture.addressMode[2] = cudaAddressModeBorder;
-	gT_par3DProjTexture.filterMode = cudaFilterModeLinear;
-	gT_par3DProjTexture.normalized = false;
-
-	cudaBindTextureToArray(gT_par3DProjTexture, array, channelDesc);
-
-	// TODO: error value?
-
-	return true;
-}
-
-
 template<unsigned int ZSIZE>
-__global__ void dev_par3D_BP(void* D_volData, unsigned int volPitch, int startAngle, int angleOffset, const SDimensions3D dims, float fOutputScale)
+__global__ void dev_par3D_BP(void* D_volData, unsigned int volPitch, cudaTextureObject_t tex, int startAngle, int angleOffset, const SDimensions3D dims, float fOutputScale)
 {
 	float* volData = (float*)D_volData;
 
@@ -125,7 +103,7 @@ __global__ void dev_par3D_BP(void* D_volData, unsigned int volPitch, int startAn
 
 			for (int idx = 0; idx < ZSIZE; ++idx) {
 
-				float fVal = tex3D(gT_par3DProjTexture, fU, fAngle, fV);
+				float fVal = tex3D<float>(tex, fU, fAngle, fV);
 				Z[idx] += fVal * fS;
 
 				fU += fCu.z;
@@ -144,7 +122,7 @@ __global__ void dev_par3D_BP(void* D_volData, unsigned int volPitch, int startAn
 }
 
 // supersampling version
-__global__ void dev_par3D_BP_SS(void* D_volData, unsigned int volPitch, int startAngle, int angleOffset, const SDimensions3D dims, int iRaysPerVoxelDim, float fOutputScale)
+__global__ void dev_par3D_BP_SS(void* D_volData, unsigned int volPitch, cudaTextureObject_t tex, int startAngle, int angleOffset, const SDimensions3D dims, int iRaysPerVoxelDim, float fOutputScale)
 {
 	float* volData = (float*)D_volData;
 
@@ -206,7 +184,7 @@ __global__ void dev_par3D_BP_SS(void* D_volData, unsigned int volPitch, int star
 				const float fU = fCu.w + fXs * fCu.x + fYs * fCu.y + fZs * fCu.z;
 				const float fV = fCv.w + fXs * fCv.x + fYs * fCv.y + fZs * fCv.z;
 
-				fVal += tex3D(gT_par3DProjTexture, fU, fAngle, fV) * fS;
+				fVal += tex3D<float>(tex, fU, fAngle, fV) * fS;
 				fZs += fSubStep;
 			}
 			fYs += fSubStep;
@@ -259,18 +237,22 @@ bool Par3DBP_Array(cudaPitchedPtr D_volumeData,
                    const SDimensions3D& dims, const SPar3DProjection* angles,
                    const SProjectorParams3D& params)
 {
-	bindProjDataTexture(D_projArray);
+	cudaTextureObject_t D_texObj;
+	if (!createTextureObject3D(D_projArray, D_texObj))
+		return false;
 
 	float fOutputScale = params.fOutputScale * params.fVolScaleX * params.fVolScaleY * params.fVolScaleZ;
 
+	bool ok = true;
+
 	for (unsigned int th = 0; th < dims.iProjAngles; th += g_MaxAngles) {
 		unsigned int angleCount = g_MaxAngles;
 		if (th + angleCount > dims.iProjAngles)
 			angleCount = dims.iProjAngles - th;
 
-		bool ok = transferConstants(angles, angleCount, params);
+		ok = transferConstants(angles, angleCount, params);
 		if (!ok)
-			return false;
+			break;
 
 		dim3 dimBlock(g_volBlockX, g_volBlockY);
 
@@ -283,23 +265,26 @@ bool Par3DBP_Array(cudaPitchedPtr D_volumeData,
 			// printf("Calling BP: %d, %dx%d, %dx%d to %p\n", i, dimBlock.x, dimBlock.y, dimGrid.x, dimGrid.y, (void*)D_volumeData.ptr); 
 			if (params.iRaysPerVoxelDim == 1) {
 				if (dims.iVolZ == 1) {
-					dev_par3D_BP<1><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), i, th, dims, fOutputScale);
+					dev_par3D_BP<1><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), D_texObj, i, th, dims, fOutputScale);
 				} else {
-					dev_par3D_BP<g_volBlockZ><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), i, th, dims, fOutputScale);
+					dev_par3D_BP<g_volBlockZ><<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), D_texObj, i, th, dims, fOutputScale);
 				}
 			} else
-				dev_par3D_BP_SS<<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), i, th, dims, params.iRaysPerVoxelDim, fOutputScale);
+				dev_par3D_BP_SS<<<dimGrid, dimBlock>>>(D_volumeData.ptr, D_volumeData.pitch/sizeof(float), D_texObj, i, th, dims, params.iRaysPerVoxelDim, fOutputScale);
 		}
 
 		// TODO: Consider not synchronizing here, if possible.
-		if (!checkCuda(cudaThreadSynchronize(), "cone_bp"))
-			return false;
+		ok = checkCuda(cudaThreadSynchronize(), "cone_bp");
+		if (!ok)
+			break;
 
 		angles = angles + angleCount;
 		// printf("%f\n", toc(t));
 
 	}
 
+	cudaDestroyTextureObject(D_texObj);
+
 	return true;
 }
 
diff --git a/cuda/3d/par3d_fp.cu b/cuda/3d/par3d_fp.cu
index 5daddc1..b2178ec 100644
--- a/cuda/3d/par3d_fp.cu
+++ b/cuda/3d/par3d_fp.cu
@@ -35,10 +35,6 @@ along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
 
 #include <cuda.h>
 
-typedef texture<float, 3, cudaReadModeElementType> texture3D;
-
-static texture3D gT_par3DVolumeTexture;
-
 namespace astraCUDA3d {
 
 static const unsigned int g_anglesPerBlock = 4;
@@ -63,24 +59,6 @@ __constant__ float gC_DetVY[g_MaxAngles];
 __constant__ float gC_DetVZ[g_MaxAngles];
 
 
-static bool bindVolumeDataTexture(const cudaArray* array)
-{
-	cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<float>();
-
-	gT_par3DVolumeTexture.addressMode[0] = cudaAddressModeBorder;
-	gT_par3DVolumeTexture.addressMode[1] = cudaAddressModeBorder;
-	gT_par3DVolumeTexture.addressMode[2] = cudaAddressModeBorder;
-	gT_par3DVolumeTexture.filterMode = cudaFilterModeLinear;
-	gT_par3DVolumeTexture.normalized = false;
-
-	cudaBindTextureToArray(gT_par3DVolumeTexture, array, channelDesc);
-
-	// TODO: error value?
-
-	return true;
-}
-
-
 // x=0, y=1, z=2
 struct DIR_X {
 	__device__ float nSlices(const SDimensions3D& dims) const { return dims.iVolX; }
@@ -89,7 +67,7 @@ struct DIR_X {
 	__device__ float c0(float x, float y, float z) const { return x; }
 	__device__ float c1(float x, float y, float z) const { return y; }
 	__device__ float c2(float x, float y, float z) const { return z; }
-	__device__ float tex(float f0, float f1, float f2) const { return tex3D(gT_par3DVolumeTexture, f0, f1, f2); }
+	__device__ float tex(cudaTextureObject_t tex, float f0, float f1, float f2) const { return tex3D<float>(tex, f0, f1, f2); }
 	__device__ float x(float f0, float f1, float f2) const { return f0; }
 	__device__ float y(float f0, float f1, float f2) const { return f1; }
 	__device__ float z(float f0, float f1, float f2) const { return f2; }
@@ -103,7 +81,7 @@ struct DIR_Y {
 	__device__ float c0(float x, float y, float z) const { return y; }
 	__device__ float c1(float x, float y, float z) const { return x; }
 	__device__ float c2(float x, float y, float z) const { return z; }
-	__device__ float tex(float f0, float f1, float f2) const { return tex3D(gT_par3DVolumeTexture, f1, f0, f2); }
+	__device__ float tex(cudaTextureObject_t tex, float f0, float f1, float f2) const { return tex3D<float>(tex, f1, f0, f2); }
 	__device__ float x(float f0, float f1, float f2) const { return f1; }
 	__device__ float y(float f0, float f1, float f2) const { return f0; }
 	__device__ float z(float f0, float f1, float f2) const { return f2; }
@@ -117,7 +95,7 @@ struct DIR_Z {
 	__device__ float c0(float x, float y, float z) const { return z; }
 	__device__ float c1(float x, float y, float z) const { return x; }
 	__device__ float c2(float x, float y, float z) const { return y; }
-	__device__ float tex(float f0, float f1, float f2) const { return tex3D(gT_par3DVolumeTexture, f1, f2, f0); }
+	__device__ float tex(cudaTextureObject_t tex, float f0, float f1, float f2) const { return tex3D<float>(tex, f1, f2, f0); }
 	__device__ float x(float f0, float f1, float f2) const { return f1; }
 	__device__ float y(float f0, float f1, float f2) const { return f2; }
 	__device__ float z(float f0, float f1, float f2) const { return f0; }
@@ -145,6 +123,7 @@ struct SCALE_NONCUBE {
 
 template<class COORD, class SCALE>
 __global__ void par3D_FP_t(float* D_projData, unsigned int projPitch,
+                           cudaTextureObject_t tex,
                            unsigned int startSlice,
                            unsigned int startAngle, unsigned int endAngle,
                            const SDimensions3D dims,
@@ -210,7 +189,7 @@ __global__ void par3D_FP_t(float* D_projData, unsigned int projPitch,
 
 		for (int s = startSlice; s < endSlice; ++s)
 		{
-			fVal += c.tex(f0, f1, f2);
+			fVal += c.tex(tex, f0, f1, f2);
 			f0 += 1.0f;
 			f1 += a1;
 			f2 += a2;
@@ -225,6 +204,7 @@ __global__ void par3D_FP_t(float* D_projData, unsigned int projPitch,
 // Supersampling version
 template<class COORD>
 __global__ void par3D_FP_SS_t(float* D_projData, unsigned int projPitch,
+                              cudaTextureObject_t tex,
                               unsigned int startSlice,
                               unsigned int startAngle, unsigned int endAngle,
                               const SDimensions3D dims, int iRaysPerDetDim,
@@ -301,7 +281,7 @@ __global__ void par3D_FP_SS_t(float* D_projData, unsigned int projPitch,
 
 		for (int s = startSlice; s < endSlice; ++s)
 		{
-			fVal += c.tex(f0, f1, f2);
+			fVal += c.tex(tex, f0, f1, f2);
 			f0 += 1.0f;
 			f1 += a1;
 			f2 += a2;
@@ -415,7 +395,9 @@ __global__ void par3D_FP_SumSqW_t(float* D_projData, unsigned int projPitch,
 
 
 bool Par3DFP_Array_internal(cudaPitchedPtr D_projData,
-                   const SDimensions3D& dims, unsigned int angleCount, const SPar3DProjection* angles,
+                   cudaTextureObject_t D_texObj,
+                   const SDimensions3D& dims,
+                   unsigned int angleCount, const SPar3DProjection* angles,
                    const SProjectorParams3D& params)
 {
 	// transfer angles to constant memory
@@ -520,29 +502,29 @@ bool Par3DFP_Array_internal(cudaPitchedPtr D_projData,
 					for (unsigned int i = 0; i < dims.iVolX; i += g_blockSlices)
 						if (params.iRaysPerDetDim == 1)
 								if (cube)
-										par3D_FP_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, scube);
+										par3D_FP_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), D_texObj, i, blockStart, blockEnd, dims, scube);
 								else
-										par3D_FP_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, snoncubeX);
+										par3D_FP_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float),  D_texObj,i, blockStart, blockEnd, dims, snoncubeX);
 						else
-							par3D_FP_SS_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeX);
+							par3D_FP_SS_t<DIR_X><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float),  D_texObj,i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeX);
 				} else if (blockDirection == 1) {
 					for (unsigned int i = 0; i < dims.iVolY; i += g_blockSlices)
 						if (params.iRaysPerDetDim == 1)
 								if (cube)
-										par3D_FP_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, scube);
+										par3D_FP_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float),  D_texObj,i, blockStart, blockEnd, dims, scube);
 								else
-										par3D_FP_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, snoncubeY);
+										par3D_FP_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float),  D_texObj,i, blockStart, blockEnd, dims, snoncubeY);
 						else
-							par3D_FP_SS_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeY);
+							par3D_FP_SS_t<DIR_Y><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float),  D_texObj,i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeY);
 				} else if (blockDirection == 2) {
 					for (unsigned int i = 0; i < dims.iVolZ; i += g_blockSlices)
 						if (params.iRaysPerDetDim == 1)
 								if (cube)
-										par3D_FP_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, scube);
+										par3D_FP_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float),  D_texObj,i, blockStart, blockEnd, dims, scube);
 								else
-										par3D_FP_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, snoncubeZ);
+										par3D_FP_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float),  D_texObj,i, blockStart, blockEnd, dims, snoncubeZ);
 						else
-							par3D_FP_SS_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float), i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeZ);
+							par3D_FP_SS_t<DIR_Z><<<dimGrid, dimBlock, 0, stream>>>((float*)D_projData.ptr, D_projData.pitch/sizeof(float),  D_texObj,i, blockStart, blockEnd, dims, params.iRaysPerDetDim, snoncubeZ);
 				}
 
 			}
@@ -569,10 +551,16 @@ bool Par3DFP(cudaPitchedPtr D_volumeData,
              const SDimensions3D& dims, const SPar3DProjection* angles,
              const SProjectorParams3D& params)
 {
+
 	// transfer volume to array
 	cudaArray* cuArray = allocateVolumeArray(dims);
 	transferVolumeToArray(D_volumeData, cuArray, dims);
-	bindVolumeDataTexture(cuArray);
+
+	cudaTextureObject_t D_texObj;
+	if (!createTextureObject3D(cuArray, D_texObj)) {
+		cudaFreeArray(cuArray);
+		return false;
+	}
 
 	bool ret;
 
@@ -584,7 +572,7 @@ bool Par3DFP(cudaPitchedPtr D_volumeData,
 		cudaPitchedPtr D_subprojData = D_projData;
 		D_subprojData.ptr = (char*)D_projData.ptr + iAngle * D_projData.pitch;
 
-		ret = Par3DFP_Array_internal(D_subprojData,
+		ret = Par3DFP_Array_internal(D_subprojData, D_texObj,
 		                             dims, iEndAngle - iAngle, angles + iAngle,
 		                             params);
 		if (!ret)
@@ -593,6 +581,8 @@ bool Par3DFP(cudaPitchedPtr D_volumeData,
 
 	cudaFreeArray(cuArray);
 
+	cudaDestroyTextureObject(D_texObj);
+
 	return ret;
 }
 
diff --git a/cuda/3d/util3d.cu b/cuda/3d/util3d.cu
index 71b5668..3dc915d 100644
--- a/cuda/3d/util3d.cu
+++ b/cuda/3d/util3d.cu
@@ -378,6 +378,28 @@ bool transferHostProjectionsToArray(const float *projData, cudaArray* array, con
 	return checkCuda(cudaMemcpy3D(&p), "transferHostProjectionsToArray 3D");
 }
 
+bool createTextureObject3D(cudaArray* array, cudaTextureObject_t& texObj)
+{
+	cudaChannelFormatDesc channelDesc =
+	    cudaCreateChannelDesc(32, 0, 0, 0, cudaChannelFormatKindFloat);
+
+	cudaResourceDesc resDesc;
+	memset(&resDesc, 0, sizeof(resDesc));
+	resDesc.resType = cudaResourceTypeArray;
+	resDesc.res.array.array = array;
+
+	cudaTextureDesc texDesc;
+	memset(&texDesc, 0, sizeof(texDesc));
+	texDesc.addressMode[0] = cudaAddressModeBorder;
+	texDesc.addressMode[1] = cudaAddressModeBorder;
+	texDesc.addressMode[2] = cudaAddressModeBorder;
+	texDesc.filterMode = cudaFilterModeLinear;
+	texDesc.readMode = cudaReadModeElementType;
+	texDesc.normalizedCoords = 0;
+
+	return checkCuda(cudaCreateTextureObject(&texObj, &resDesc, &texDesc, NULL), "createTextureObject3D");
+}
+
 
 
 float dotProduct3D(cudaPitchedPtr data, unsigned int x, unsigned int y,