Merge branch 'adv-raycaster' into 'development'

Adv raycaster

See merge request !108

application/campvisapplication.cpp

@@ -235,14 +235,8 @@ namespace campvis {
 
             _mainWindow->deinit();
             QuadRenderer::deinit();
-
-            // deinit OpenGL and cgt
-            cgt::deinitGL();
         }
 
-        // MainWindow dtor needs a valid CampVisApplication, so we need to call it here instead of during destruction.
-        delete _mainWindow;
-
         // now delete everything in the right order:
         for (std::vector<PipelineRecord>::iterator it = _pipelines.begin(); it != _pipelines.end(); ++it) {
             delete it->_painter;
@@ -252,6 +246,12 @@ namespace campvis {
             delete *it;
         }
 
+        {
+            // Deinit everything OpenGL using the local context.
+            cgt::GLContextScopedLock lock(_localContext);
+            cgt::deinitGL();
+        }
+
         GLJobProc.stop();
         OpenGLJobProcessor::deinit();
         SimpleJobProcessor::deinit();

core/classification/abstracttransferfunction.h

@@ -90,6 +90,13 @@ namespace campvis {
          */
         virtual size_t getDimensionality() const = 0;
 
+        /**
+         * Returns the intensity domain where this TF has it's non-transparent parts.
+         * I.e. the minimum and the maximum intensity being opaque.
+         * \return  The intensity domain where this TF has it's non-transparent parts.
+         */
+        virtual cgt::vec2 getVisibilityDomain() const = 0;
+
         /**
          * Binds the transfer function OpenGL texture to the given texture and sets up uniforms.
          * \note    Calling thread must have a valid OpenGL context.

core/classification/genericgeometrytransferfunction.h

@@ -59,6 +59,13 @@ namespace campvis {
          * Destructor, make sure to delete the OpenGL texture beforehand by calling deinit() with a valid OpenGL context!
          */
         virtual ~GenericGeometryTransferFunction();
+
+        /**
+         * Returns the intensity domain where this TF has it's non-transparent parts.
+         * I.e. the minimum and the maximum intensity being opaque.
+         * \return  cgt::vec2(0.f, 1.f)
+         */
+        virtual cgt::vec2 getVisibilityDomain() const;
         
         /**
          * Initializes the Shader, hence, this methods has to be called from a thread with a valid OpenGL context!
@@ -128,6 +135,20 @@ namespace campvis {
     campvis::GenericGeometryTransferFunction<T>::~GenericGeometryTransferFunction() {
     }
 
+    template<class T>
+    cgt::vec2 campvis::GenericGeometryTransferFunction<T>::getVisibilityDomain() const {
+        if (_geometries.empty())
+            return cgt::vec2(-1.f, -1.f);
+        else {
+            cgt::vec2 minmax(1.f, 0.f);
+            for (size_t i = 0; i < _geometries.size(); ++i) {
+                minmax.x = std::min(minmax.x, _geometries[i]->getIntensityDomain().x);
+                minmax.y = std::max(minmax.y, _geometries[i]->getIntensityDomain().y);
+            }
+            return minmax;
+        }
+    }
+
     template<class T>
     void campvis::GenericGeometryTransferFunction<T>::initShader() {
         _shader = ShdrMgr.load("core/glsl/passthrough.vert", "core/glsl/passthrough.frag", "");

core/classification/simpletransferfunction.cpp

@@ -103,5 +103,9 @@ namespace campvis {
         return _rightColor;
     }
 
+    cgt::vec2 SimpleTransferFunction::getVisibilityDomain() const {
+        return cgt::vec2(0.f, 1.f);
+    }
+
 
 }
\ No newline at end of file

core/classification/simpletransferfunction.h

@@ -58,6 +58,13 @@ namespace campvis {
          */
         virtual size_t getDimensionality() const;
 
+        /**
+         * Returns the intensity domain where this TF has it's non-transparent parts.
+         * I.e. the minimum and the maximum intensity being opaque.
+         * \return  cgt::vec2(0.f, 1.f)
+         */
+        virtual cgt::vec2 getVisibilityDomain() const;
+
         const cgt::col4& getLeftColor() const;
         void setLeftColor(const cgt::col4& color);
         const cgt::col4& getRightColor() const;

core/classification/tfgeometry1d.cpp

@@ -157,4 +157,8 @@ namespace campvis {
         _keyPoints.insert(lb, kp);
     }
 
+    cgt::vec2 TFGeometry1D::getIntensityDomain() const {
+        return cgt::vec2(_keyPoints.front()._position, _keyPoints.back()._position);
+    }
+
 }
\ No newline at end of file

core/classification/tfgeometry1d.h

@@ -110,6 +110,11 @@ namespace campvis {
          */
         void addKeyPoint(float position, const cgt::col4& color);
 
+        /**
+         * Returns the intensity domain of this TFGeometry1D.
+         */
+        cgt::vec2 getIntensityDomain() const;
+
         /**
          * Creates a simple quad geometry for the given interval.
          * A quad geometry consists of two KeyPoints.

core/classification/tfgeometry2d.cpp

@@ -129,6 +129,10 @@ namespace campvis {
         std::sort(_keyPoints.begin(), _keyPoints.end(), KeyPointSorter(_center._position));
     }
 
+    cgt::vec2 TFGeometry2D::getIntensityDomain() const {
+        return cgt::vec2(_keyPoints.front()._position.x, _keyPoints.back()._position.x);
+    }
+
 
 
 }
\ No newline at end of file

core/classification/tfgeometry2d.h

@@ -75,6 +75,11 @@ namespace campvis {
          * \return 
          */
         std::vector<KeyPoint>& getKeyPoints();
+        
+        /**
+         * Returns the intensity domain of this TFGeometry1D.
+         */
+        cgt::vec2 getIntensityDomain() const;
 
         /**
          * Renders this transfer function geometry to the current active OpenGL context.

core/glsl/tools/texture3d.frag

@@ -59,8 +59,9 @@ vec4 textureToWorld(in TextureParameters3D texParams, in vec4 texCoords) {
  * \param   texCoords   texture coordinates
  * \return  \a texCoords transformes to woorld coordinates.
  */
-vec4 textureToWorld(in TextureParameters3D texParams, in vec3 texCoords) {
-    return textureToWorld(texParams, vec4(texCoords, 1.0));
+vec3 textureToWorld(in TextureParameters3D texParams, in vec3 texCoords) {
+    vec4 v = textureToWorld(texParams, vec4(texCoords, 1.0));
+    return v.xyz;
 }
 
 /**

modules/preprocessing/processors/glimageresampler.cpp

@@ -44,11 +44,13 @@ namespace campvis {
         , p_inputImage("InputImage", "Input Image", "", DataNameProperty::READ)
         , p_outputImage("OutputImage", "Output Image", "GlImageResampler.out", DataNameProperty::WRITE)
         , p_resampleScale("ResampleScale", "Resampling Scale", .5f, .01f, 10.f)
+        , p_targetSize("TargetSize", "Size of Resampled Image", cgt::ivec3(128), cgt::ivec3(1), cgt::ivec3(1024))
         , _shader(0)
     {
-        addProperty(p_inputImage);
+        addProperty(p_inputImage, INVALID_RESULT | INVALID_PROPERTIES);
         addProperty(p_outputImage);
-        addProperty(p_resampleScale);
+        addProperty(p_resampleScale, INVALID_RESULT | INVALID_PROPERTIES);
+        addProperty(p_targetSize);
     }
 
     GlImageResampler::~GlImageResampler() {
@@ -73,7 +75,7 @@ namespace campvis {
 
         if (img != 0) {
             cgt::vec3 originalSize(img->getSize());
-            cgt::ivec3 resampledSize(cgt::ceil(originalSize * p_resampleScale.getValue()));
+            const cgt::ivec3& resampledSize = p_targetSize.getValue();
 
             cgt::TextureUnit inputUnit;
             inputUnit.activate();
@@ -115,4 +117,14 @@ namespace campvis {
         }
     }
 
+    void GlImageResampler::updateProperties(DataContainer& dataContainer) {
+        ImageRepresentationGL::ScopedRepresentation img(dataContainer, p_inputImage.getValue());
+
+        if (img != 0) {
+            p_targetSize.setMaxValue(cgt::ivec3(img->getSize()) * int(p_resampleScale.getMaxValue()));
+            p_targetSize.setValue(cgt::ivec3(cgt::vec3(img->getSize()) * p_resampleScale.getValue()));
+        }
+    }
+
 }
+

modules/preprocessing/processors/glimageresampler.h

@@ -72,11 +72,13 @@ namespace campvis {
         DataNameProperty p_outputImage;     ///< ID for output gradient volume
 
         FloatProperty p_resampleScale;      ///< Resampling Scale
+        IVec3Property p_targetSize;         ///< Size of resampled image
 
     protected:
         /// \see AbstractProcessor::updateResult
         virtual void updateResult(DataContainer& dataContainer);
-
+        virtual void updateProperties(DataContainer& dataContainer);
+        
         cgt::Shader* _shader;               ///< Shader for resampling
 
         static const std::string loggerCat_;
@@ -85,3 +87,4 @@ namespace campvis {
 }
 
 #endif // GLIMAGERESAMPLER_H__
+

modules/vis/glsl/advoptimizedraycaster.frag

@@ -33,6 +33,8 @@ layout(location = 2) out vec4 out_FHN;       ///< outgoing fragment first hit no
 #include "tools/texture3d.frag"
 #include "tools/transferfunction.frag"
 
+#include "modules/vis/glsl/voxelhierarchy.frag"
+
 uniform vec2 _viewportSizeRCP;
 uniform float _jitterStepSizeMultiplier;
 
@@ -47,7 +49,7 @@ uniform sampler2D _exitPoints;
 uniform sampler2D _exitPointsDepth;
 uniform TextureParameters2D _exitParams;
 
-// DRR volume
+// Input volume
 uniform sampler3D _volume;
 uniform TextureParameters3D _volumeTextureParams;
 
@@ -55,183 +57,35 @@ uniform TextureParameters3D _volumeTextureParams;
 uniform sampler1D _transferFunction;
 uniform TFParameters1D _transferFunctionParams;
 
-
-// BBV Lookup volume
-uniform usampler2D _vvTexture;
-uniform int _vvVoxelSize;
-uniform int _vvVoxelDepth;
-uniform int _vvMaxMipMapLevel;
-
-
 uniform LightSource _lightSource;
 uniform vec3 _cameraPosition;
 
 uniform float _samplingStepSize;
-
-#ifdef INTERSECTION_REFINEMENT
-bool _inVoid = false;
-#endif
-
-#ifdef ENABLE_SHADOWING
-uniform float _shadowIntensity;
-#endif
-
-// TODO: copy+paste from Voreen - eliminate or improve.
 const float SAMPLING_BASE_INTERVAL_RCP = 200.0;
 
-// converting the sample coordinates [0, 1]^3 to integer voxel coordinate in l-level of the voxelized texture. 
-ivec3 voxelToBrick(in vec3 voxel, int level) {
-    int res = 1 << level;
-
-    return ivec3(floor(voxel.x / (_vvVoxelSize * res)), floor(voxel.y / (_vvVoxelSize * res)), floor(voxel.z / _vvVoxelDepth));
-}
-
-// Looks up the l mipmap level of voxelized texture of the render data and returns the value in samplePosition coordinates.
-// samplePosition is in texture coordiantes [0, 1]
-// level is the level in the hierarchy
-uint lookupInBbv(in vec3 samplePosition, in int level) {
-    ivec3 byte = voxelToBrick(samplePosition * _volumeTextureParams._size, level);
-    uint texel = uint(texelFetch(_vvTexture, byte.xy, level).r);
-
-    return texel;
-}
-
-// Converts the ray to bitmask which can be intersect with voxelized volume
-void rayMapToMask(in vec3 samplePos, in vec3 entryPoint, in vec3 direction, in int level, out uint mask, inout float tFar) {
-    int res = 1 << level ;							//< the number of voxels from voxelized volume data which are mapped into 
-                                                    //< voxelized volume in l-th level of hierarchy 
-
-    // calculating the voxel which the ray's origin will start from.
-    // BoxLlf is the lower-left corner of the voxel.
-    // BoxUrb is the upper-right corner of the voxel.
-    // NOTE: the z-direction coordinate is not important. So, we simply put BoxLlf.z = 0 and BoxUrb.z = 1. Because we are 
-    //		 going to intersect the raybitmask with the bitmask of the voxel in z-direction to find the intersection point.
-    vec3 brickVoxel = floor((samplePos * _volumeTextureParams._size) / vec3(_vvVoxelSize * res, _vvVoxelSize * res, _vvVoxelDepth)) * vec3(_vvVoxelSize * res, _vvVoxelSize * res, _vvVoxelDepth);
-    vec3 boxLlf = vec3((brickVoxel * _volumeTextureParams._sizeRCP).xy, 0.0f);
-    vec3 boxUrb = vec3((boxLlf + (_volumeTextureParams._sizeRCP * vec3(_vvVoxelSize * res, _vvVoxelSize * res, _vvVoxelDepth))).xy, 1.0f);
-
-    // here, tmin and tmax for the ray which is intersecting with the voxel is computed. 
-    vec3 tMin = (boxLlf - entryPoint) / direction;
-    vec3 tMax = (boxUrb - entryPoint) / direction;
-
-    // As tmin and tmax values are not necessary tmin and tmax values and value should be resorted, to find the 
-    // exit point of the voxel we just compute minimum value of maximum values in tmin and tmax.
-    vec3 tFarVec = max(tMin, tMax);
-    tFar = min(tFar, min(tFarVec.x, min(tFarVec.y, tFarVec.z)));
-
-    // Here we are going to calculate the entry and exit point from the voxel by calculated tFar.
-    int bit1 = clamp(int(floor(samplePos.z * 32)), 0, 31);
-    int bit2 = clamp(int(floor((entryPoint.z + tFar * direction.z) * 32)), 0, 31);
-
-    // setting the bits of the bitmask values.
-    mask = bitfieldInsert(uint(0x0), uint(0xffffffff), max((min(bit1, bit2)), 0), min(abs(bit2 - bit1) + 1, 32));
-    
-    // HACK HACK HACK
-    // right now, the bitmask generation for ray is not working correctly.
-    mask = 0xffffffff;
-}
-
-// intersects the ray bitmask and voxelized data bitmap.
-// bitray: ray bitmask,
-// texel: the texel coordinate in the voxelized volume data,
-// level: the level of the texture,
-// intersectionBitmask: the result bitmask of intersecting the data in texel of the voxelized volume data and ray bitmask.
-// it returns wherther there was an intersection or not.
-bool intersectBits(in uint bitRay, in vec3 texel, in int level, out uint intersectionBitmask) {
-   // Fetch bitmask from hierarchy and compute intersection via bitwise AND
-   intersectionBitmask = (bitRay & lookupInBbv(texel, level));
-   return (intersectionBitmask != uint(0));
-}
-
-// 
-bool intersectHierarchy(in int level, in vec3 samplePos, in vec3 entryPoint, in vec3 direction, inout float tFar, out uint intersectionBitmask) {
-   uint rayBitmask = uint(0);
-
-   // Mapping the ray to a bitmask
-   rayMapToMask(samplePos, entryPoint, direction, level, rayBitmask, tFar);
-
-   // Intersect the ray bitmask with current pixel
-   return intersectBits(rayBitmask, samplePos, level, intersectionBitmask);
-}
-
 /**
  * Performs the raycasting and returns the final fragment color.
  */
 vec4 performRaycasting(in vec3 entryPoint, in vec3 exitPoint, in vec2 texCoords) {
     vec4 result = vec4(0.0);
     vec3 direction = exitPoint - entryPoint;
-    float tNear = 0.0;
-    float tFar = length(direction);
-    direction = normalize(direction);
+    
+    // Adjust direction a bit to prevent division by zero
+    direction.x = (abs(direction.x) < 0.0000001) ? 0.0000001 : direction.x;
+    direction.y = (abs(direction.y) < 0.0000001) ? 0.0000001 : direction.y;
+    direction.z = (abs(direction.z) < 0.0000001) ? 0.0000001 : direction.z;
+
+    OFFSET = (0.25 / (1 << _vhMaxMipMapLevel)); //< offset value used to avoid self-intersection or previous voxel intersection.
 
     jitterEntryPoint(entryPoint, direction, _samplingStepSize * _jitterStepSizeMultiplier);
 
     float firstHitT = -1.0f;
-  
- 
-    // compute sample position
-    vec3 samplePosition = entryPoint.rgb + tNear * direction;
-
-    bool intersectionFound = false;
-    uint intersectionBitmask = uint(0);
-    int level = min(3, _vvMaxMipMapLevel);			//< start from a level of hierarchy. We are trying to start from one of the levels in middle (not coarse and not best).
-    float offset = (0.0025 / int(1 << int(_vvMaxMipMapLevel))) / length(direction);	//< offset value used to avoid self-intersection or previous voxel intersection.
-    int i = 0;
-    float newTfar = 1.0f;							//< tFar calculated for the point where the ray is going out of current texel in voxelized data.
-    while(!intersectionFound && (tNear <= tFar) && (level >= 0) && (level < _vvMaxMipMapLevel) && (i < 100)) {
-        newTfar = tFar;								//< because it will be compared to calculated tFar in intersectHierarchy function.
-
-        // intersects the ray with the hierarchy of voxelized volume texture.
-        if(intersectHierarchy(level, samplePosition, entryPoint, direction, newTfar, intersectionBitmask)) {
-            // if the level is 0, the intersection is found, otherwise check the intersection in a higher level of hierarchy.
-            intersectionFound = (level == 0);               
-            level--;
-        } else {
-            // The beginning part of the ray is cut and going to a coarser level of hierarchy.
-            tNear = newTfar + offset;
-            samplePosition = entryPoint + tNear * direction;
-            level++;                
-        }
-
-        i++;
-    }
-
-              
-    if(!intersectionFound) {
-        result = vec4(0, 0, 0, 0);
-        return result;
-    }
-
-    bool tFarIntersectionFound = false;
-    uint tFarIntersectionBitmask = uint(0);
-    int tFarLevel = min(3, _vvMaxMipMapLevel);
-    int k = 0;
-    float tFarDec = 0.0;
-    float newTnear = 1.0f;
-    vec3 tFarSamplePosition = exitPoint - tFarDec * direction;
-    while(!tFarIntersectionFound && (tFarDec <= tFar) && (tFarLevel >= 0) && (tFarLevel < _vvMaxMipMapLevel) && (k < 100)) {
-        newTnear = tFar;
-
-        if(intersectHierarchy(tFarLevel, tFarSamplePosition, exitPoint, -direction, newTnear, tFarIntersectionBitmask)) {
-            tFarIntersectionFound = (tFarLevel == 0);               
-            tFarLevel--;
-        } else {
-            tFarDec = newTnear + offset;
-            tFarSamplePosition = exitPoint - tFarDec * direction;
-            tFarLevel++;
 
-                
-        }
+    float tNear = clipFirstHitpoint(entryPoint, direction, 0.0, 1.0);
+    float tFar = 1.0 - clipFirstHitpoint(exitPoint, -direction, 0.0, 1.0);
 
-        k++;
-    }
-              
-    if(!tFarIntersectionFound) {
-        result = vec4(0, 0, 0, 0);
-        return result;
-    }
-
-    tFar -= tFarDec;
+    // compute sample position
+    vec3 samplePosition = entryPoint.rgb + tNear * direction;
 
     while (tNear < tFar) {
         vec3 samplePosition = entryPoint.rgb + tNear * direction;
@@ -240,46 +94,13 @@ vec4 performRaycasting(in vec3 entryPoint, in vec3 exitPoint, in vec2 texCoords)
         float intensity = texture(_volume, samplePosition).r;
         vec4 color = lookupTF(_transferFunction, _transferFunctionParams, intensity);
 
-#ifdef INTERSECTION_REFINEMENT
-        if (color.a <= 0.0) {
-            // we're within void, make the steps bigger
-            _inVoid = true;
-        } else {
-            if (_inVoid) {
-                float formerT = t - _samplingStepSize;
-
-                // we just left the void, perform intersection refinement
-                for (float stepPower = 0.5; stepPower > 0.1; stepPower /= 2.0) {
-                    // compute refined sample position
-                    float newT = formerT + _samplingStepSize * stepPower;
-                    vec3 newSamplePosition = entryPoint.rgb + newT * direction;
-
-                    // lookup refined intensity + TF
-                    float newIntensity = texture(_volume, newSamplePosition).r;
-                    vec4 newColor = lookupTF(_transferFunction, _transferFunctionParams, newIntensity);
-
-                    if (newColor.a <= 0.0) {
-                        // we're back in the void - look on the right-hand side
-                        formerT = newT;
-                    }
-                    else {
-                        // we're still in the matter - look on the left-hand side
-                        samplePosition = newSamplePosition;
-                        color = newColor;
-                        t -= _samplingStepSize * stepPower;
-                    }
-                }
-                _inVoid = false;
-            }
-        }
-#endif
-
         // perform compositing
         if (color.a > 0.0) {
 #ifdef ENABLE_SHADING
             // compute gradient (needed for shading and normals)
             vec3 gradient = computeGradient(_volume, _volumeTextureParams, samplePosition);
-            color.rgb = calculatePhongShading(textureToWorld(_volumeTextureParams, samplePosition).xyz, _lightSource, _cameraPosition, gradient, color.rgb, color.rgb, vec3(1.0, 1.0, 1.0));
+            vec4 worldPos = _volumeTextureParams._textureToWorldMatrix * vec4(samplePosition, 1.0); // calling textureToWorld here crashes Intel HD driver and nVidia driver in debug mode, hence, let's calc it manually...
+            color.rgb = calculatePhongShading(worldPos.xyz / worldPos.w, _lightSource, _cameraPosition, gradient, color.rgb);
 #endif
 
             // accomodate for variable sampling rates
@@ -299,7 +120,7 @@ vec4 performRaycasting(in vec3 entryPoint, in vec3 exitPoint, in vec2 texCoords)
         if (result.a > 0.975) {
             result.a = 1.0;
             tNear = tFar;
-
+        
                 
         }
 

modules/vis/glsl/contextpreservingraycaster.frag

@@ -103,8 +103,8 @@ vec4 performRaycasting(in vec3 entryPoint, in vec3 exitPoint, in vec2 texCoords)
             // compute gradient (needed for shading and normals)
             vec3 gradient = computeGradient(_volume, _volumeTextureParams, samplePosition) * 10;
 
-            color.rgb = calculatePhongShading(textureToWorld(_volumeTextureParams, samplePosition).xyz, _lightSource, _cameraPosition, gradient, color.rgb, color.rgb, vec3(1.0, 1.0, 1.0));
-            float SI = getPhongShadingIntensity(textureToWorld(_volumeTextureParams, samplePosition).xyz, _lightSource, _cameraPosition, gradient);
+            color.rgb = calculatePhongShading(textureToWorld(_volumeTextureParams, samplePosition), _lightSource, _cameraPosition, gradient, color.rgb, color.rgb, vec3(1.0, 1.0, 1.0));
+            float SI = getPhongShadingIntensity(textureToWorld(_volumeTextureParams, samplePosition), _lightSource, _cameraPosition, gradient);
 
             // accomodate for variable sampling rates
             color.a = 1.0 - pow(1.0 - color.a, _samplingStepSize * SAMPLING_BASE_INTERVAL_RCP);

modules/vis/glsl/hierarchyrenderer.frag

@@ -25,7 +25,7 @@
 in vec3 ex_TexCoord;        ///< incoming texture coordinate
 in vec4 ex_Position;        ///< incoming texture coordinate
 
-out uint result;
+out uvec4 result;
  
 #include "tools/texture2d.frag"
 #include "tools/texture3d.frag"
@@ -42,47 +42,56 @@ uniform TFParameters1D _transferFunctionParams;
 uniform float _inverseTexSizeX; // 1.0 / mipmapLevelResolution
 uniform float _inverseTexSizeY; // 1.0 / mipmapLevelResolution
 
-uniform uint _voxelSize;
-uniform uint _voxelDepth;
+uniform int _brickSize;
+uniform int _brickDepth;
+uniform vec3 _hierarchySize;
+uniform vec2 _tfDomain;
 
 void main() {
-    result = uint(0);
+    result = uvec4(0);
 
-    vec3 startIdx = vec3(0, 0, 0);
-    vec3 endIdx = _volumeTextureParams._size - vec3(1, 1, 1);
+    // For each element in the uvec4 bitmask:
+    for (int e = 0; e < 4; ++e) {
+        // For each bit of the uint value, the corresponding area in the volume is checked, whether it's transparent or not. The bits are set accordingly.
+        for (int d = 0; d < 32; d++) {
+            vec3 llf = vec3(gl_FragCoord.xy / _hierarchySize.xy, float((e * 32) + d + 0.5) / 128.0);
+            //ivec3 llf = ivec3((gl_FragCoord.x - 0.5) * _brickSize, (gl_FragCoord.y - 0.5) * _brickSize, );
+            bool hasData = false;
 
-    // For each bit of the int value, the area is checked to find voxels or not. If some voxels were found, the bit is set or cleared.
-    for(int d = 0; d < 32; d++){
+            /** 
+             * For each bit _brickSize number of voxels in x and y direction and _brickDepth number of voxel in z-direction should be considered.
+             * Also, to make sure that all voxels are considered, one offset voxel is considered in x, y, and z boundary from each side so that
+             * each side will consider 2 more voxels.
+             */
+            for (int z = -1; z < _brickDepth + 1; ++z) {
+                for (int y = -1; y < _brickSize + 1; ++y) {
+                    for (int x = -1; x < _brickSize + 1; ++x) {
+                        vec3 addendum = (vec3(x, y, z) / _volumeTextureParams._size);
+                        vec3 texCoord = clamp(llf + addendum, 0.0, 1.0);
+                        //float intensity = mapIntensityToTFDomain(_transferFunctionParams._intensityDomain, texture(_volume, texCoord).r);
+                        ivec3 voxel = ivec3(texCoord * _volumeTextureParams._size);
+                        float intensity = mapIntensityToTFDomain(_transferFunctionParams._intensityDomain, texelFetch(_volume, voxel, 0).r);
+                        //float intensity = texture(_volume, voxel).r;
+                        //vec4 color = lookupTF(_transferFunction, _transferFunctionParams, intensity);
+                        //if (color.a > 0) {
+                        //}
 
-        
-        vec3 samplePosition = vec3((gl_FragCoord.x - 0.5) * _voxelSize, (gl_FragCoord.y - 0.5) * _voxelSize, d * _voxelDepth);
-
-        bool hasData = false;
-        float intensity;
-        vec4 color;
-
-        /** For each bit _voxelSize number of voxels in x and y direction and _voxelDepth number of voxel in z-direction should be considered.
-         * Also, to make sure that all voxels are considered, one offset voxel is considered in x, y, and z boundary from each side so that
-         * each side will consider 2 more voxels.
-         */
-        for(int i = 0; i < (_voxelSize + 2) * (_voxelSize + 2) * (_voxelDepth + 2); i++){
-
-            // the offset value calculates which voxel should be checked.
-            ivec3 offset = ivec3(i % (_voxelSize + 2) - 1, (i / (_voxelSize + 2)) % (_voxelSize + 2) - 1, i / ((_voxelSize + 2) * (_voxelSize + 2)) - 1);
-            
-            vec3 s = clamp(samplePosition + vec3(0.5, 0.5, 0.5) + offset, startIdx, endIdx);
-            intensity = texture(_volume, s * _volumeTextureParams._sizeRCP).r;
-            color = lookupTF(_transferFunction, _transferFunctionParams, intensity);
-
-            // if there was any data in the volume data in that voxel, set the bit.
-            if(color.a > 0){
-                hasData = true;
-                break;
+                        // if there was any data in the volume data in that voxel, set the bit.
+                        if (intensity >= _tfDomain.x && intensity <= _tfDomain.y) {