<div dir="ltr">No worries...so we just checked the math the cell to point offset matrix seems to be right but the delta seems wrong. Can you change this <div><br></div><div>> delta[0] = this->Extents[1] - this->Extents[0];<br>> delta[1] = this->Extents[3] - this->Extents[2];<br>> delta[2] = this->Extents[5] - this->Extents[4];<br><br></div><div>to this</div><div><br></div><div>delta[0] = this->Extents[1] - this->Extents[0] + 1;<br>delta[1] = this->Extents[3] - this->Extents[2] + 1;<br>delta[2] = this->Extents[5] - this->Extents[4] + 1;</div><div><br></div><div>Thanks,</div><div><br><br></div><div><br></div><div><br><div><div><br><div class="gmail_quote"><div dir="ltr">On Thu, Oct 5, 2017 at 10:48 AM Elvis Stansvik <<a href="mailto:elvis.stansvik@orexplore.com" target="_blank">elvis.stansvik@orexplore.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">2017-10-05 16:34 GMT+02:00 Aashish Chaudhary <<a href="mailto:aashish.chaudhary@kitware.com" target="_blank">aashish.chaudhary@kitware.com</a>>:<br>
> Thanks we are looking into it, earlier we moved the math from the shader to<br>
> a matrix that now we are forming in C++ so at this point we are looking at<br>
> the equation carefully. But yes, I thought you may have found a bug.<br>
<br>
Thanks for having a look. I suspect my std::max(..) is just a<br>
workaround and not a proper fix.<br>
<br>
Me and a colleague had another look, but I think we were both too<br>
tired to figure out if something was fundamentally wrong with the<br>
equations. Too late in the day :)<br>
<br>
Elvis<br>
<br>
><br>
> - aashish<br>
><br>
> On Thu, Oct 5, 2017 at 8:39 AM Elvis Stansvik <<a href="mailto:elvis.stansvik@orexplore.com" target="_blank">elvis.stansvik@orexplore.com</a>><br>
> wrote:<br>
>><br>
>> 2017-10-05 10:38 GMT+02:00 Elvis Stansvik <<a href="mailto:elvis.stansvik@orexplore.com" target="_blank">elvis.stansvik@orexplore.com</a>>:<br>
>> > 2017-10-04 18:12 GMT+02:00 Aashish Chaudhary<br>
>> > <<a href="mailto:aashish.chaudhary@kitware.com" target="_blank">aashish.chaudhary@kitware.com</a>>:<br>
>> >> Got it, thanks!<br>
>> ><br>
>> > Although I'm not familiar at all with the code, I suspect this is what<br>
>> > happens with my 200x200x2 volume:<br>
>> ><br>
>> > In void<br>
>> > vtkOpenGLGPUVolumeRayCastMapper::vtkInternal::ComputeCellToPointMatrix()<br>
>> > {<br>
>> > this->CellToPointMatrix->Identity();<br>
>> > this->AdjustedTexMin[0] = this->AdjustedTexMin[1] =<br>
>> > this->AdjustedTexMin[2] = 0.0f;<br>
>> > this->AdjustedTexMin[3] = 1.0f;<br>
>> > this->AdjustedTexMax[0] = this->AdjustedTexMax[1] =<br>
>> > this->AdjustedTexMax[2] = 1.0f;<br>
>> > this->AdjustedTexMax[3] = 1.0f;<br>
>> ><br>
>> > if (!this->Parent->CellFlag) // point data<br>
>> > {<br>
>> > float delta[3];<br>
>> > delta[0] = this->Extents[1] - this->Extents[0];<br>
>> > delta[1] = this->Extents[3] - this->Extents[2];<br>
>> > delta[2] = this->Extents[5] - this->Extents[4];<br>
>> ><br>
>> > float min[3];<br>
>> > min[0] = 0.5f / delta[0];<br>
>> > min[1] = 0.5f / delta[1];<br>
>> > min[2] = 0.5f / delta[2];<br>
>> ><br>
>> > float range[3]; // max - min<br>
>> > range[0] = (delta[0] - 0.5f) / delta[0] - min[0];<br>
>> > range[1] = (delta[1] - 0.5f) / delta[1] - min[1];<br>
>> > range[2] = (delta[2] - 0.5f) / delta[2] - min[2];<br>
>> ><br>
>> > this->CellToPointMatrix->SetElement(0, 0, range[0]); // Scale diag<br>
>> > this->CellToPointMatrix->SetElement(1, 1, range[1]);<br>
>> > this->CellToPointMatrix->SetElement(2, 2, range[2]);<br>
>> > this->CellToPointMatrix->SetElement(0, 3, min[0]); // t vector<br>
>> > this->CellToPointMatrix->SetElement(1, 3, min[1]);<br>
>> > this->CellToPointMatrix->SetElement(2, 3, min[2]);<br>
>> ><br>
>> > // Adjust limit coordinates for texture access.<br>
>> > float const zeros[4] = {0.0f, 0.0f, 0.0f, 1.0f}; // GL tex min<br>
>> > float const ones[4] = {1.0f, 1.0f, 1.0f, 1.0f}; // GL tex max<br>
>> > this->CellToPointMatrix->MultiplyPoint(zeros, this->AdjustedTexMin);<br>
>> > this->CellToPointMatrix->MultiplyPoint(ones, this->AdjustedTexMax);<br>
>> > }<br>
>> > }<br>
>> ><br>
>> > delta[2] = 1 ==> min[2] = 0.5 ==> range[2] = 0<br>
>> ><br>
>> > So the scale factor in the third dimension will be 0 in<br>
>> > CellToPointMatrix, right? I guess this could be what accounts for the<br>
>> > strange rendering?<br>
>> ><br>
>> > Could/should this code be updated to take into account when the extent<br>
>> > delta is just 1, so that volumes 2 data points thick can be rendered?<br>
>><br>
>> I've only tested it briefly, but this seems to give the correct<br>
>> rendering even with a 2 point thick input:<br>
>><br>
>> diff --git a/Rendering/VolumeOpenGL2/vtkOpenGLGPUVolumeRayCastMapper.cxx<br>
>> b/Rendering/VolumeOpenGL2/vtkOpenGLGPUVolumeRayCastMapper.cxx<br>
>> index 8767bd0..88bbcba 100644<br>
>> --- a/Rendering/VolumeOpenGL2/vtkOpenGLGPUVolumeRayCastMapper.cxx<br>
>> +++ b/Rendering/VolumeOpenGL2/vtkOpenGLGPUVolumeRayCastMapper.cxx<br>
>> @@ -1495,9 +1495,9 @@ void<br>
>> vtkOpenGLGPUVolumeRayCastMapper::vtkInternal::ComputeCellToPointMatrix()<br>
>> if (!this->Parent->CellFlag) // point data<br>
>> {<br>
>> float delta[3];<br>
>> - delta[0] = this->Extents[1] - this->Extents[0];<br>
>> - delta[1] = this->Extents[3] - this->Extents[2];<br>
>> - delta[2] = this->Extents[5] - this->Extents[4];<br>
>> + delta[0] = std::max(2, this->Extents[1] - this->Extents[0]);<br>
>> + delta[1] = std::max(2, this->Extents[3] - this->Extents[2]);<br>
>> + delta[2] = std::max(2, this->Extents[5] - this->Extents[4]);<br>
>><br>
>> float min[3];<br>
>> min[0] = 0.5f / delta[0];<br>
>><br>
>> This ensures the texture does not collapse into 0 when the extent<br>
>> delta is 1. Would you consider this a correct fix?<br>
>><br>
>> Elvis<br>
>><br>
>> ><br>
>> > Elvis<br>
>> ><br>
>> >><br>
>> >> On Wed, Oct 4, 2017 at 12:07 PM Elvis Stansvik<br>
>> >> <<a href="mailto:elvis.stansvik@orexplore.com" target="_blank">elvis.stansvik@orexplore.com</a>> wrote:<br>
>> >>><br>
>> >>> 2017-10-04 17:43 GMT+02:00 Aashish Chaudhary<br>
>> >>> <<a href="mailto:aashish.chaudhary@kitware.com" target="_blank">aashish.chaudhary@kitware.com</a>>:<br>
>> >>> > thanks..I just wanted to make sure that we remove that as one of the<br>
>> >>> > culprit. I am wondering if it is something to do with the geometry<br>
>> >>> > (bbox) or<br>
>> >>> > something to do with some edge case (calculation that is going off<br>
>> >>> > because<br>
>> >>> > of very small delta). Although in the past we have tried very think<br>
>> >>> > volumes<br>
>> >>> ><br>
>> >>> > Is it possible for you to send a sample dataset to us?<br>
>> >>><br>
>> >>> My test case includes the data used. It's simply a 200x200x2 volume<br>
>> >>> filled with all data points set to 0.01.<br>
>> >>><br>
>> >>> Elvis<br>
>> >>><br>
>> >>> ><br>
>> >>> > On Wed, Oct 4, 2017 at 10:09 AM Elvis Stansvik<br>
>> >>> > <<a href="mailto:elvis.stansvik@orexplore.com" target="_blank">elvis.stansvik@orexplore.com</a>> wrote:<br>
>> >>> >><br>
>> >>> >> Hi Aashish,<br>
>> >>> >><br>
>> >>> >> 2017-10-04 15:47 GMT+02:00 Aashish Chaudhary<br>
>> >>> >> <<a href="mailto:aashish.chaudhary@kitware.com" target="_blank">aashish.chaudhary@kitware.com</a>>:<br>
>> >>> >> > Please make sure that your sampling distance is less than < 2.<br>
>> >>> >> > You<br>
>> >>> >> > can<br>
>> >>> >> > set<br>
>> >>> >> > the sampling distance manually in the API.<br>
>> >>> >><br>
>> >>> >> Hm. I don't think it's related to sampling distance. In the test<br>
>> >>> >> case<br>
>> >>> >> I'm not setting sampling distance at all, so I think it should<br>
>> >>> >> default<br>
>> >>> >> to 1, and then be automatically adjusted during interaction.<br>
>> >>> >><br>
>> >>> >> In any case, even if I amend the test case with e.g.:<br>
>> >>> >><br>
>> >>> >> volumeMapper->AutoAdjustSampleDistancesOff();<br>
>> >>> >> volumeMapper->SetSampleDistance(0.2);<br>
>> >>> >> ...<br>
>> >>> >> volumeMapper2->AutoAdjustSampleDistancesOff();<br>
>> >>> >> volumeMapper2->SetSampleDistance(0.2);<br>
>> >>> >><br>
>> >>> >> The result is the same.<br>
>> >>> >><br>
>> >>> >> Elvis<br>
>> >>> >><br>
>> >>> >> ><br>
>> >>> >> > On Wed, Oct 4, 2017 at 8:23 AM Elvis Stansvik<br>
>> >>> >> > <<a href="mailto:elvis.stansvik@orexplore.com" target="_blank">elvis.stansvik@orexplore.com</a>><br>
>> >>> >> > wrote:<br>
>> >>> >> >><br>
>> >>> >> >> Hi all,<br>
>> >>> >> >><br>
>> >>> >> >> In the test case below, I render one 200x200x2 volume and one<br>
>> >>> >> >> 200x200x3 volume, both with the same settings.<br>
>> >>> >> >><br>
>> >>> >> >> The 200x200x3 volume looks as expected (the light gray one in<br>
>> >>> >> >> the<br>
>> >>> >> >> attached screenshot), while the 200x200x2 one looks strange (too<br>
>> >>> >> >> dark,<br>
>> >>> >> >> and with some kind of gradient artifact).<br>
>> >>> >> >><br>
>> >>> >> >> Is there some limitation on rendering a volume that is only 2<br>
>> >>> >> >> thick<br>
>> >>> >> >> in<br>
>> >>> >> >> one of the dimensions? I can understand why a volume of<br>
>> >>> >> >> thickness 1<br>
>> >>> >> >> can't be rendered properly, since VTK must have values to<br>
>> >>> >> >> interpolate<br>
>> >>> >> >> between, but 2 should be OK right?<br>
>> >>> >> >><br>
>> >>> >> >> Thanks for any help.<br>
>> >>> >> >><br>
>> >>> >> >> Elvis<br>
>> >>> >> >><br>
>> >>> >> >><br>
>> >>> >> >> main.cpp:<br>
>> >>> >> >><br>
>> >>> >> >> #include <algorithm><br>
>> >>> >> >><br>
>> >>> >> >> #include <vtkColorTransferFunction.h><br>
>> >>> >> >> #include <vtkGenericOpenGLRenderWindow.h><br>
>> >>> >> >> #include <vtkGPUVolumeRayCastMapper.h><br>
>> >>> >> >> #include <vtkImageData.h><br>
>> >>> >> >> #include <vtkNew.h><br>
>> >>> >> >> #include <vtkPiecewiseFunction.h><br>
>> >>> >> >> #include <vtkProperty.h><br>
>> >>> >> >> #include <vtkRenderer.h><br>
>> >>> >> >> #include <vtkRenderWindow.h><br>
>> >>> >> >> #include <vtkRenderWindowInteractor.h><br>
>> >>> >> >> #include <vtkVolume.h><br>
>> >>> >> >> #include <vtkVolumeProperty.h><br>
>> >>> >> >><br>
>> >>> >> >> int main(int argc, char *argv[])<br>
>> >>> >> >> {<br>
>> >>> >> >> vtkNew<vtkColorTransferFunction> colorFunction;<br>
>> >>> >> >> colorFunction->AddRGBPoint(0.0, 0.0, 0.0, 0.0);<br>
>> >>> >> >> colorFunction->AddRGBPoint(1.0, 0.0, 0.0, 0.0);<br>
>> >>> >> >><br>
>> >>> >> >> vtkNew<vtkPiecewiseFunction> opacityFunction;<br>
>> >>> >> >> opacityFunction->AddPoint(0.0, 0.0);<br>
>> >>> >> >> opacityFunction->AddPoint(1.0, 1.0);<br>
>> >>> >> >><br>
>> >>> >> >> // A 200x200x2 volume<br>
>> >>> >> >> vtkNew<vtkImageData> imageData;<br>
>> >>> >> >> imageData->SetExtent(0, 199, 0, 199, 0, 1);<br>
>> >>> >> >> imageData->AllocateScalars(VTK_FLOAT, 1);<br>
>> >>> >> >> std::fill_n(static_cast<float<br>
>> >>> >> >> *>(imageData->GetScalarPointer()),<br>
>> >>> >> >> 200 * 200 * 2, 0.01);<br>
>> >>> >> >><br>
>> >>> >> >> vtkNew<vtkGPUVolumeRayCastMapper> volumeMapper;<br>
>> >>> >> >> volumeMapper->SetInputData(imageData.Get());<br>
>> >>> >> >><br>
>> >>> >> >> vtkNew<vtkVolumeProperty> volumeProperty;<br>
>> >>> >> >> volumeProperty->SetScalarOpacity(opacityFunction.Get());<br>
>> >>> >> >> volumeProperty->SetColor(colorFunction.Get());<br>
>> >>> >> >> volumeProperty->ShadeOff();<br>
>> >>> >> >><br>
>> >>> >> >> vtkNew<vtkVolume> volume;<br>
>> >>> >> >> volume->SetMapper(volumeMapper.Get());<br>
>> >>> >> >> volume->SetProperty(volumeProperty.Get());<br>
>> >>> >> >> volume->SetPosition(120, 0, 0);<br>
>> >>> >> >><br>
>> >>> >> >> // A 200x200x3 volume<br>
>> >>> >> >> vtkNew<vtkImageData> imageData2;<br>
>> >>> >> >> imageData2->SetExtent(0, 199, 0, 199, 0, 2);<br>
>> >>> >> >> imageData2->AllocateScalars(VTK_FLOAT, 1);<br>
>> >>> >> >> std::fill_n(static_cast<float<br>
>> >>> >> >> *>(imageData2->GetScalarPointer()),<br>
>> >>> >> >> 200 * 200 * 3, 0.01);<br>
>> >>> >> >><br>
>> >>> >> >> vtkNew<vtkGPUVolumeRayCastMapper> volumeMapper2;<br>
>> >>> >> >> volumeMapper2->SetInputData(imageData2.Get());<br>
>> >>> >> >><br>
>> >>> >> >> vtkNew<vtkVolumeProperty> volumeProperty2;<br>
>> >>> >> >> volumeProperty2->SetScalarOpacity(opacityFunction.Get());<br>
>> >>> >> >> volumeProperty2->SetColor(colorFunction.Get());<br>
>> >>> >> >> volumeProperty2->ShadeOff();<br>
>> >>> >> >><br>
>> >>> >> >> vtkNew<vtkVolume> volume2;<br>
>> >>> >> >> volume2->SetMapper(volumeMapper2.Get());<br>
>> >>> >> >> volume2->SetProperty(volumeProperty2.Get());<br>
>> >>> >> >> volume2->SetPosition(-120, 0, 0);<br>
>> >>> >> >><br>
>> >>> >> >> vtkNew<vtkRenderer> renderer;<br>
>> >>> >> >> renderer->AddVolume(volume.Get());<br>
>> >>> >> >> renderer->AddVolume(volume2.Get());<br>
>> >>> >> >> renderer->SetBackground(1.0, 1.0, 1.0);<br>
>> >>> >> >><br>
>> >>> >> >> // Render with "plain" render window / interactor<br>
>> >>> >> >> vtkNew<vtkRenderWindow> window;<br>
>> >>> >> >> window->SetMultiSamples(0);<br>
>> >>> >> >> window->AddRenderer(renderer.Get());<br>
>> >>> >> >><br>
>> >>> >> >> vtkNew<vtkRenderWindowInteractor> interactor;<br>
>> >>> >> >> interactor->SetRenderWindow(window.Get());<br>
>> >>> >> >> interactor->Start();<br>
>> >>> >> >><br>
>> >>> >> >> return 0;<br>
>> >>> >> >> }<br>
>> >>> >> >><br>
>> >>> >> >><br>
>> >>> >> >> CMakeLists.txt:<br>
>> >>> >> >><br>
>> >>> >> >> cmake_minimum_required(VERSION 3.1)<br>
>> >>> >> >><br>
>> >>> >> >> project(TestCase)<br>
>> >>> >> >><br>
>> >>> >> >> find_package(VTK 8.0 COMPONENTS<br>
>> >>> >> >> vtkCommonCore<br>
>> >>> >> >> vtkCommonDataModel<br>
>> >>> >> >> vtkCommonExecutionModel<br>
>> >>> >> >> vtkCommonMath<br>
>> >>> >> >> vtkFiltersSources<br>
>> >>> >> >> vtkGUISupportQt<br>
>> >>> >> >> vtkInteractionStyle<br>
>> >>> >> >> vtkRenderingCore<br>
>> >>> >> >> vtkRenderingOpenGL2<br>
>> >>> >> >> vtkRenderingVolume<br>
>> >>> >> >> vtkRenderingVolumeOpenGL2<br>
>> >>> >> >> REQUIRED<br>
>> >>> >> >> )<br>
>> >>> >> >><br>
>> >>> >> >> add_executable(TestCase main.cpp)<br>
>> >>> >> >><br>
>> >>> >> >> target_link_libraries(TestCase PUBLIC<br>
>> >>> >> >> vtkCommonCore<br>
>> >>> >> >> vtkCommonDataModel<br>
>> >>> >> >> vtkCommonExecutionModel<br>
>> >>> >> >> vtkCommonMath<br>
>> >>> >> >> vtkFiltersSources<br>
>> >>> >> >> vtkInteractionStyle<br>
>> >>> >> >> vtkRenderingCore<br>
>> >>> >> >> vtkRenderingOpenGL2<br>
>> >>> >> >> vtkRenderingVolume<br>
>> >>> >> >> vtkRenderingVolumeOpenGL2<br>
>> >>> >> >> )<br>
>> >>> >> >><br>
>> >>> >> >> target_include_directories(TestCase PUBLIC<br>
>> >>> >> >> ${VTK_INCLUDE_DIRS}<br>
>> >>> >> >> )<br>
>> >>> >> >><br>
>> >>> >> >> target_compile_definitions(TestCase PUBLIC<br>
>> >>> >> >> ${VTK_DEFINITIONS}<br>
>> >>> >> >> )<br>
>> >>> >> >><br>
>> >>> >> >> set_target_properties(TestCase PROPERTIES<br>
>> >>> >> >> CXX_STANDARD 14<br>
>> >>> >> >> CXX_STANDARD_REQUIRED ON<br>
>> >>> >> >> )<br>
>> >>> >> >> _______________________________________________<br>
>> >>> >> >> Powered by <a href="http://www.kitware.com" rel="noreferrer" target="_blank">www.kitware.com</a><br>
>> >>> >> >><br>
>> >>> >> >> Visit other Kitware open-source projects at<br>
>> >>> >> >> <a href="http://www.kitware.com/opensource/opensource.html" rel="noreferrer" target="_blank">http://www.kitware.com/opensource/opensource.html</a><br>
>> >>> >> >><br>
>> >>> >> >> Please keep messages on-topic and check the VTK FAQ at:<br>
>> >>> >> >> <a href="http://www.vtk.org/Wiki/VTK_FAQ" rel="noreferrer" target="_blank">http://www.vtk.org/Wiki/VTK_FAQ</a><br>
>> >>> >> >><br>
>> >>> >> >> Search the list archives at:<br>
>> >>> >> >> <a href="http://markmail.org/search/?q=vtkusers" rel="noreferrer" target="_blank">http://markmail.org/search/?q=vtkusers</a><br>
>> >>> >> >><br>
>> >>> >> >> Follow this link to subscribe/unsubscribe:<br>
>> >>> >> >> <a href="http://public.kitware.com/mailman/listinfo/vtkusers" rel="noreferrer" target="_blank">http://public.kitware.com/mailman/listinfo/vtkusers</a><br>
</blockquote></div></div></div></div></div>