[Insight-developers] Get/SetPixel() with taking an offsetparameter
Miller, James V (GE, Research)
millerjv at crd.ge.com
Wed Dec 20 16:33:34 EST 2006
Gaetan,
You might also want to look at using std::deque instead of std::vector. I probably overuse std::vector and should use std::deque more. Deque doesn't use contiguous memory, so when it needs to increase its size it just allocates another block (it is basically a link list of fixed length vectors).
Using a deque, you won't have to prescan so you can allocate vectors of the exact length. You can be lazy and use push_back() and have simpler code. Plus you don't pay the log2-type allocation penalty of a vector as the deque grows. So you should get a performance increase over not preallocating the vectors.
So you might be able to use an std::map<PixelType, std::deque<IndexType> > or an std::map<PixelType, std::deque<OffsetValueType> >.
If we do introduce an image type with non-continguous memory, I will push that we maintain Image with an API that provides a one-dimensional view of ComputeIndex()/ComputeOffset(). What I would want to "remove" is the access to GetBufferPointer() and hence prevent people writing code equivalent to
x = img->GetBufferPointer() + img->ComputeOffset(index);
*x = val;
But still allow a 1D access
offset = img->ComputeOffset(index);
img->SetPixel(offset, val);
Or
for (i=0; i < numPixel; ++i)
img->SetPixel(I, val);
So from a programming perspective, you can think of Image as having a continguous block of memory. But internally, we'd have the freedom to divide the image up differently.
We can do this if we modify the iterators and the IO. We'd have to pay a double dereference penalty to access a pixel. And we'd upset anyone out in the world that is relying on GetBufferPointer(). But we may be able to do this using new image types and traits for creating iterators, etc....
Those are my thoughts.
Jim
-----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: Wednesday, December 20, 2006 4:37 AM
To: Miller, James V (GE, Research); Luis Ibanez
Cc: insight-developers at itk.org
Subject: Re: [Insight-developers] Get/SetPixel() with taking an offsetparameter
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=s
>>> ummary),
>>> ...) 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
_______________________________________________
Insight-developers mailing list
Insight-developers at itk.org
http://www.itk.org/mailman/listinfo/insight-developers
More information about the Insight-developers
mailing list