<div dir="ltr">Hi Bradley, thank for your response. I have already check the itk module you send me, in fact I use it as a reference to create my own, yours and mine basically do the same. I will try to use yours in the same 3D image and see if there is any problem similar to mine.<div>
<br></div><div>Thanks a lot.</div></div><div class="gmail_extra"><br><br><div class="gmail_quote">On Mon, Aug 19, 2013 at 8:31 PM, Bradley Lowekamp <span dir="ltr"><<a href="mailto:blowekamp@mail.nih.gov" target="_blank">blowekamp@mail.nih.gov</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Hello,<br>
<br>
I am not entirely sure what your filter is trying to do, nor what's the problem you are encountering.<br>
<br>
I have an itk module which has a filter to compute the GLCM for a neighborhood around each pixel. This may be what you are looking for, or at least give you some ideas on how to use the components.<br>
<br>
As it has been a while since I last used this filter, and there is a lack of rigorous testing, I'd recommend writing a test to verify the results for your usage before relaying on the output.<br>
<br>
Good luck,<br>
Brad<br>
<br>
[1] <a href="https://github.com/blowekamp/itkTextureAnalysis/blob/master/include/itkTextureFeatureImageFilter.h" target="_blank">https://github.com/blowekamp/itkTextureAnalysis/blob/master/include/itkTextureFeatureImageFilter.h</a><br>
<div><div class="h5"><br>
On Aug 19, 2013, at 8:46 PM, Fabian Torres <<a href="mailto:dae.wong@gmail.com">dae.wong@gmail.com</a>> wrote:<br>
<br>
> Hi all. I am trying to implement a filter that gives local texture analysis. For this I implement a itkImagetoImageFilter that calculates coocurrence texture descriptor for each pixel using itkScalarImageToCoocurrenceMatrizFilter.<br>
><br>
> I know this is a slow implementation, but the real problem I have is that when I try to use this with a 3D image, even a small one of 43x35x61, the computer runs out of memory. I do not know if I my implementation is wrong. But it seems that for each pixel I´m using a lot of memory.<br>
><br>
> I hope someone can help me with this. Here I leave the code for the GenerateData() function.<br>
><br>
> Thanks.<br>
><br>
> this->AllocateOutputs();<br>
><br>
> ProgressReporter progress( this, 0,<br>
> this->GetInput()->GetRequestedRegion().GetNumberOfPixels(), 100 );<br>
><br>
> typedef MirrorPadImageFilter<InputImageType,InputImageType> PadFilterType;<br>
> typename PadFilterType::Pointer padFilter = PadFilterType::New();<br>
><br>
> int boundSize = floor(this->m_RegionSize[1]/2);<br>
><br>
> RegionSizeType bound;<br>
> bound.Fill(boundSize);<br>
><br>
> padFilter->SetPadBound(bound);<br>
> padFilter->SetInput(this->GetInput());<br>
> padFilter->Update();<br>
><br>
> typename InputImageType::Pointer padImage = padFilter->GetOutput();<br>
><br>
> typedef Statistics::ScalarImageToCooccurrenceMatrixFilter<InputImageType><br>
> GLCMGeneratorType;<br>
> typename GLCMGeneratorType::Pointer glcmGenerator = GLCMGeneratorType::New();<br>
><br>
> glcmGenerator->SetNumberOfBinsPerAxis(this->m_NumberofBins);<br>
> glcmGenerator->SetPixelValueMinMax(this->m_PixelMinValue,m_PixelMaxValue);<br>
><br>
> typedef typename InputImageType::OffsetType OffsetType;<br>
><br>
> OffsetType offset1;<br>
> offset1[0] = 1;<br>
> offset1[1] = 0;<br>
><br>
> OffsetType offset2;<br>
> offset2[0] = 1;<br>
> offset2[1] = -1;<br>
><br>
> OffsetType offset3;<br>
> offset3[0] = 0;<br>
> offset3[1] = -1;<br>
><br>
> OffsetType offset4;<br>
> offset4[0] = -1;<br>
> offset4[1] = -1;<br>
><br>
> typedef typename GLCMGeneratorType::OffsetVector OffsetVectorType;<br>
> typename OffsetVectorType::Pointer offsets = OffsetVectorType::New();<br>
><br>
> offsets->reserve(4);<br>
> offsets->InsertElement(0,offset1);<br>
> offsets->InsertElement(1,offset2);<br>
> offsets->InsertElement(2,offset3);<br>
> offsets->InsertElement(3,offset4);<br>
><br>
> glcmGenerator->SetOffsets(offsets);<br>
><br>
> typedef typename GLCMGeneratorType::HistogramType HistogramType;<br>
> typedef Statistics::HistogramToTextureFeaturesFilter<HistogramType><br>
> TextureFeaturesType;<br>
> typename TextureFeaturesType::Pointer textureFeatures;<br>
><br>
> if(this->m_TextureFeature<=8 && this->m_TextureFeature>=0){<br>
> textureFeatures = TextureFeaturesType::New();<br>
> }else{<br>
> textureFeatures = NULL;<br>
> }<br>
><br>
> typedef ExtractImageFilter<InputImageType,InputImageType> ExtractFilterType;<br>
> typename ExtractFilterType::Pointer extractFilter = ExtractFilterType::New();<br>
><br>
> typename InputImageType::RegionType region;<br>
> typename InputImageType::IndexType regionIndex;<br>
><br>
> region.SetSize(this->m_RegionSize);<br>
><br>
> extractFilter->SetInput(padImage);<br>
><br>
> ImageRegionIteratorType itPad(padImage, padImage->GetLargestPossibleRegion());<br>
> itPad.GoToBegin();<br>
><br>
> ImageRegionIteratorType itOut(this->GetOutput(),<br>
> this->GetOutput()->GetLargestPossibleRegion());<br>
> itOut.GoToBegin();<br>
><br>
> while(!itPad.IsAtEnd()){<br>
><br>
> regionIndex = itPad.GetIndex();<br>
> region.SetIndex(regionIndex);<br>
><br>
> if(padImage->GetLargestPossibleRegion().IsInside(region)){<br>
><br>
> extractFilter->SetExtractionRegion(region);<br>
> extractFilter->UpdateLargestPossibleRegion();<br>
><br>
> glcmGenerator->SetInput(extractFilter->GetOutput());<br>
> glcmGenerator->UpdateLargestPossibleRegion();<br>
><br>
> if(this->m_TextureFeature == 0){<br>
> textureFeatures->SetInput(glcmGenerator->GetOutput());<br>
> textureFeatures->UpdateLargestPossibleRegion();<br>
> this->GetOutput()->SetPixel(itOut.GetIndex(),<br>
> textureFeatures->GetFeature(TextureFeaturesType::Energy));<br>
> }else if(this->m_TextureFeature == 1){<br>
> textureFeatures->SetInput(glcmGenerator->GetOutput());<br>
> textureFeatures->UpdateLargestPossibleRegion();<br>
> this->GetOutput()->SetPixel(itOut.GetIndex(),<br>
> textureFeatures->GetFeature(TextureFeaturesType::Entropy));<br>
> }else if(this->m_TextureFeature == 2){<br>
> textureFeatures->SetInput(glcmGenerator->GetOutput());<br>
> textureFeatures->UpdateLargestPossibleRegion();<br>
> this->GetOutput()->SetPixel(itOut.GetIndex(),<br>
> textureFeatures->GetFeature(TextureFeaturesType::Correlation));<br>
> }else if(this->m_TextureFeature == 3){<br>
> textureFeatures->SetInput(glcmGenerator->GetOutput());<br>
> textureFeatures->UpdateLargestPossibleRegion();<br>
> this->GetOutput()->SetPixel(itOut.GetIndex(),<br>
> textureFeatures->GetFeature(TextureFeaturesType::InverseDifferenceMoment));<br>
> }else if(this->m_TextureFeature == 4){<br>
> textureFeatures->SetInput(glcmGenerator->GetOutput());<br>
> textureFeatures->UpdateLargestPossibleRegion();<br>
> this->GetOutput()->SetPixel(itOut.GetIndex(),<br>
> textureFeatures->GetFeature(TextureFeaturesType::Inertia));<br>
> }else if(this->m_TextureFeature == 5){<br>
> textureFeatures->SetInput(glcmGenerator->GetOutput());<br>
> textureFeatures->UpdateLargestPossibleRegion();<br>
> this->GetOutput()->SetPixel(itOut.GetIndex(),<br>
> textureFeatures->GetFeature(TextureFeaturesType::ClusterShade));<br>
> }else if(this->m_TextureFeature == 6){<br>
> textureFeatures->SetInput(glcmGenerator->GetOutput());<br>
> textureFeatures->UpdateLargestPossibleRegion();<br>
> this->GetOutput()->SetPixel(itOut.GetIndex(),<br>
> textureFeatures->GetFeature(TextureFeaturesType::ClusterProminence));<br>
> }else if(this->m_TextureFeature == 7){<br>
> textureFeatures->SetInput(glcmGenerator->GetOutput());<br>
> textureFeatures->UpdateLargestPossibleRegion();<br>
> this->GetOutput()->SetPixel(itOut.GetIndex(),<br>
> textureFeatures->GetFeature(TextureFeaturesType::HaralickCorrelation));<br>
> }else if(this->m_TextureFeature == 8){<br>
><br>
> typename HistogramType::ConstPointer hist = glcmGenerator->GetOutput();<br>
><br>
> typename HistogramType::ConstIterator it = hist->Begin();<br>
><br>
> float mean = 0;<br>
> while(it != hist->End()){<br>
> mean += it.GetFrequency();<br>
> ++it;<br>
> }<br>
><br>
> float nBins = hist->GetSize(0);<br>
> nBins *= nBins;<br>
><br>
> mean /= nBins;<br>
><br>
> float variance = 0;<br>
> it = hist->Begin();<br>
> while( it != hist->End()){<br>
> variance += pow(it.GetFrequency()-mean,2);<br>
> ++it;<br>
> }<br>
><br>
> variance /= (nBins - 1);<br>
><br>
> this->GetOutput()->SetPixel(itOut.GetIndex(),variance);<br>
><br>
> hist = NULL;<br>
><br>
> }<br>
><br>
> ++itOut;<br>
><br>
> progress.CompletedPixel();<br>
><br>
> }<br>
><br>
> ++itPad;<br>
</div></div>> _____________________________________<br>
> Powered by <a href="http://www.kitware.com" 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" target="_blank">http://www.kitware.com/opensource/opensource.html</a><br>
><br>
> Kitware offers ITK Training Courses, for more information visit:<br>
> <a href="http://www.kitware.com/products/protraining.php" target="_blank">http://www.kitware.com/products/protraining.php</a><br>
><br>
> Please keep messages on-topic and check the ITK FAQ at:<br>
> <a href="http://www.itk.org/Wiki/ITK_FAQ" target="_blank">http://www.itk.org/Wiki/ITK_FAQ</a><br>
><br>
> Follow this link to subscribe/unsubscribe:<br>
> <a href="http://www.itk.org/mailman/listinfo/insight-users" target="_blank">http://www.itk.org/mailman/listinfo/insight-users</a><br>
<br>
</blockquote></div><br><br clear="all"><div><br></div>-- <br>Fabián Torres Robles<br>Maestria en Ciencias en Ingeniería Electrónica<br>Ingeniería en Sistemas Electrónicos<br>tel. 58081280, 0445534661338<br>e-mail <a href="mailto:fabian.trobles@gmail.com" target="_blank">fabian.trobles@gmail.com</a>, <a href="mailto:dae.wong@gmail.com" target="_blank">dae.wong@gmail.com</a> <div style="padding:0px;margin-left:0px;margin-top:0px;overflow:hidden;word-wrap:break-word;color:black;font-size:10px;text-align:left;line-height:130%">
</div>
</div>