[Insight-developers] Get/SetPixel() with taking an offsetparameter
Gaetan Lehmann
gaetan.lehmann at jouy.inra.fr
Wed Dec 20 04:36:51 EST 2006
Here is what I have tried about std::map< PixelType, std::vector<
IndexType > >:
(a) - make a first pass to compute the number of pixels for all the values
- create the map, and reserve the exact size which will be required
for all the vector
- fill the vectors
(b) replace std::map< PixelType, std::vector< IndexType > > by std::map<
PixelType, std::vector< OffsetValueType > >, the conversion being done by
the ComputeOffset() and computeIndex() methods.
both (a) and (b) where done to decrease the memory usage - and it works
very well regarding that.
More surprising, both methods also give better performance. With both, my
test run in 8.2 sec instead of 10.7 before the changes.
I think that's because the vectors no more need to increase their size,
and that it avoid the cost of this increase, and because manipulating 3
times less memory compensate the cost of the conversion.
The only problem is the much compex code, with additional ComputeOffset()
and ComputeIndex() methods in it because of (b), and the 2 more steps done
in (a).
And there is still the cost of the conversion index <-> offset that I
would like to avoid
Gaetan
On Tue, 19 Dec 2006 22:48:52 +0100, Miller, James V (GE, Research)
<millerjv at crd.ge.com> wrote:
> Wrt to the std::map< PixelType, std::vector< IndexType > >, I have
> thought about this a few times and always verred away because of the
> memory requirements.
>
> If we are adopting these high memory implementations, we may want to
> have an option to switch to a lower memory albeit slower version.
>
> I just had a thought that maybe some of these algorithms that use a map
> as above could be done in several iterations. Less than the number of
> pixel intensity values but more than one. Say you divide the intensity
> space into k bands, say you build the map with the top 10% of the
> pixels, then do whatever processing you can do, then build the map with
> the next 10%, etc.
>
>
> -----Original Message-----
> From: insight-developers-bounces+millerjv=crd.ge.com at itk.org
> [mailto:insight-developers-bounces+millerjv=crd.ge.com at itk.org] On
> Behalf Of Gaetan Lehmann
> Sent: Tuesday, December 12, 2006 6:13 AM
> To: Luis Ibanez
> Cc: insight-developers at itk.org
> Subject: Re: [Insight-developers] Get/SetPixel() with taking an
> offsetparameter
>
>
> Hi,
>
>>
>> Hi Gaetan,
>>
>> That's an interesting option.
>>
>> However, please note that a better way of approaching this problem is
>> to create an Iterator that is specific to the way you want to walk in
>> your image.
>>
>> It seems that in your algorithm you are already defining a clever way
>> of finding the offsets that will lead you to the next pixel of
>> interest. What we probably want to explore is the possibility of
>> creating a new Image Iterator that will do the equivalent job.
>>
>
> What is done in ReconstructionByDilationImageFilter currently work only
> with shape iterators, so having the same operators on the other
> iterators would be great. I'm not sure that a dedicated iterator is
> needed.
>
>> In general you want to avoid having direct access to the image buffer.
>> One of the reasons for it is that soon we will be introducing in ITK
>> an alternative memory layout where the pixel data is not stored in a
>> contiguous block of memory. Instead, the pixel data may be stored as a
>> collection of disconnected slices.
>>
>> Having a GetPixel()/SetPixel() methods with direct access to the
>> buffer data will not work in this memory layout. We already will have
>> to manage the modification of the existing ComputOffset() method.
>>
>
> I'm not sure that modifying the Get/SetPixel() method is the best way to
> go, but what is proposed has no extra cost, because the ComputeIndex()
> and
> ComputeOffset() are already exposed, and there is nothing more done that
> using these methods.
>
>> You will also have a harder time making sure that the algorithm works
>> in N-Dimensional images. It is in general, easier, to delegate to
>> Iterators the complex details of dealing with the N-Dimensionality of
>> the image.
>
> It will not be harder, because an iterator is used to get the offsets
> :-) Offsets are only a cheaper representation of Index, and the
> conversion can be easily done with the ComputeIndex() and
> ComputeOffset() methods.
>
>>
>> What is the specific way in which you need to walk through the image
>> in the morphological reconstruction filters that you are preparing ?
>>
>
> A common step - at least common enough to make me implement it several
> times - is to sort the pixels by there value to process the highest value
> first for example. The data structure used to do that is:
>
> std::map< PixelType, std::vector< IndexType > >
>
> after have filled this data structure with all the pixels in the image,
> it's easy to process all the pixels in the image in the right order.
> If we take care of the memory allocation of the vectors (which require an
> extra iteration step), for a big image (1024 x 1024 x 80), this structure
> take about 1024*1024*80*4*3/1024^2 = 960 MB on a 32 bits system - twice
> more a 64bits one. When using offsets, the cost is reduced by the image
> dimension, 3 in that case = 320 MB. The difference is quite significant
> :-)
>
> This example is "only" an optimization problem, as the pixels can be
> processed in the right order without that, by iterating over the image as
> many time as we have pixel values in the image (and one more time to find
> the maximum).
>
> The next case is more difficult.
> The component tree representation (also called min tree or max tree)
> don't
> store the image in a array, but in a tree of connected components in the
> image. All the nodes of that tree are associated to a list of pixels in
> the image - a list of itk::Index. That's no more an algorithm
> optimization, but a data representation problem, and the itk way of
> designing a pixel location make it highly inefficient compared to what
> can
> be done in other image analysis libs.
>
> It is already possible to use offsets instead of itk::Index. That's what
> I'm doing in the component trees, for memory efficiency (please, don't
> make the ComputeOffset() and ComputeIndex() method deprecated !).
> However,
> the permanent conversion itk::Index <-> offset seem to be an important
> extra cost which should be avoided.
> I have not measured the extra cost of the conversion, so perhaps I'm
> wrong, and the conversion cost is not significant compared to the rest of
> the algorithms.
> Given the comments about performance in the code, I think that the
> measures of the execution time of those conversion have already been
> done.
> Can you give a pointer to those results ?
>
> Thanks,
>
> Gaetan
>
>
>>
>> Please let us know,
>>
>>
>> Thanks
>>
>>
>> Luis
>>
>>
>>
>> ----------------------
>> Gaetan Lehmann wrote:
>>> Hi,
>>> Currently, the GetPixel() and SetPixel() methods of the itk::Image
>>> class are taking a itk::Index parameter to know which pixel is
>>> concerned in the image.
>>> The itk::Image also have the ComputeIndex() and ComputeOffset() methods
>>> to manipulate directly the offset in the buffer, and to convert this
>>> offset from/to an index.
>>> Manipulating indexes in an algorithm (like in the current
>>> morphological reconstruction filters, in the component tree
>>> contribution I'm preparing
>>> (http://voxel.jouy.inra.fr/darcsweb/darcsweb.cgi?r=componentTree;a=summary),
>>> ...) would be a lot more efficient, both in memory and in execution
>>> time, by using the offsets instead of the indexes: store an itk::Index
>>> cost 3 long (with dim=3), while an offset is only 1, and there is no
>>> need to convert the offset to an index and the index to an offset.
>>> Unfortunately, there is no Get/SetPixel() method able to work with an
>>> offset.
>>> I think that it would be a good idea to implement them as follow (not
>>> tested yet, that's an example):
>>> in itkImage.h:
>>> void SetPixel(const OffsetValueType &offset, const TPixel& value)
>>> {
>>> (*m_Buffer)[offset] = value;
>>> }
>>> const TPixel& GetPixel(const OffsetValueType &offset) const
>>> {
>>> return ( (*m_Buffer)[offset] );
>>> }
>>> TPixel& GetPixel(const OffsetValueType &offset)
>>> {
>>> return ( (*m_Buffer)[offset] );
>>> }
>>> TPixel & operator[](const OffsetValueType &offset)
>>> { return this->GetPixel(offset); }
>>> const TPixel& operator[](const OffsetValueType &offset) const
>>> { return this->GetPixel(offset); }
>>> void SetPixel(const IndexType &index, const TPixel& value)
>>> {
>>> typename Superclass::OffsetValueType offset =
>>> this->ComputeOffset(index);
>>> this->SetPixel( offset, value );
>>> }
>>> const TPixel& GetPixel(const IndexType &index) const
>>> {
>>> typename Superclass::OffsetValueType offset =
>>> this->ComputeOffset(index);
>>> return this->GetPixel( offset );
>>> }
>>> TPixel& GetPixel(const IndexType &index)
>>> {
>>> typename Superclass::OffsetValueType offset =
>>> this->ComputeOffset(index);
>>> return this->GetPixel( offset );
>>> }
>>> TPixel & operator[](const IndexType &index)
>>> { return this->GetPixel(index); }
>>> const TPixel& operator[](const IndexType &index) const
>>> { return this->GetPixel(index); }
>>> in itkImageAdaptor.h:
>>> void SetPixel(const IndexType &index, const PixelType & value)
>>> { m_PixelAccessor.Set( m_Image->GetPixel(index), value ); }
>>> PixelType GetPixel(const IndexType &index) const
>>> { return m_PixelAccessor.Get( m_Image->GetPixel(index) ); }
>>> PixelType operator[](const IndexType &index) const
>>> { return m_PixelAccessor.Get( m_Image->GetPixel(index) ); }
>>> void SetPixel(const OffsetValueType &offset, const PixelType & value)
>>> { m_PixelAccessor.Set( m_Image->GetPixel(offset), value ); }
>>> PixelType GetPixel(const OffsetValueType &offset) const
>>> { return m_PixelAccessor.Get( m_Image->GetPixel(offset) ); }
>>> PixelType operator[](const OffsetValueType &offset) const
>>> { return m_PixelAccessor.Get( m_Image->GetPixel(offset) ); }
>>> Does it seem reasonable to do that ?
>>> Thanks,
>>> Gaetan
>>>
>
>
>
--
Gaëtan Lehmann
Biologie du Développement et de la Reproduction
INRA de Jouy-en-Josas (France)
tel: +33 1 34 65 29 66 fax: 01 34 65 29 09
http://voxel.jouy.inra.fr
More information about the Insight-developers
mailing list