[Insight-users] map the colours from the eigenvector

[Seth Gilchrist]

Hi Alba,
A program I wrote which I think should do something like you want is below.
I'm working in Linux, but it should compile for you.

#include <iostream>
#include <stdlib.h>
#include "itkVector.h"
#include "itkImage.h"
#include "itkImageFileWriter.h"
#include "itkRGBPixel.h"
#include "itkImageRegionIterator.h"
#include <time.h>
#include <limits>

int main(int argc, char** argv)
{
    if ( argc < 2 || argc > 2){ //Check for correct number of inputs
        std::cout<<"Usage:"<<std::endl;
        std::cout<<argv[0]<<" [Output File Name] "<<std::endl<<std::endl;
        return EXIT_FAILURE;
    }
    /* Create the Images.*/
    typedef double    VectorComponentType;
    typedef unsigned short    RGBComponentType; // many image formats have
limited pixel types
    typedef itk::Vector< VectorComponentType, 3>    VectorPixelType;
    typedef itk::RGBPixel< RGBComponentType >        RGBPixelType;

    typedef itk::Image< VectorPixelType, 3 >        VectorImageType;
    typedef itk::Image< RGBPixelType, 3 >            RGBImageType;

    VectorImageType::RegionType                region;
    VectorImageType::RegionType::SizeType    regionSize;
    VectorImageType::RegionType::IndexType  regionStart;
    VectorImageType::PointType                origin;
    VectorImageType::SpacingType            spacing;

    regionSize[0] = 100;
    regionSize[1] = 100;
    regionSize[2] = 100;
    regionStart[0] = 0;
    regionStart[1] = 0;
    regionStart[2] = 0;
    region.SetSize(regionSize);
    region.SetIndex(regionStart);
    origin[0] = 0;
    origin[1] = 0;
    origin[2] = 0;
    spacing[0] = 1;
    spacing[1] = 1;
    spacing[2] = 1;

    VectorImageType::Pointer    vectorIm = VectorImageType::New();
    RGBImageType::Pointer        rgbIm = RGBImageType::New();

    vectorIm->SetRegions(region);
    vectorIm->SetOrigin(origin);
    vectorIm->SetSpacing(spacing);
    rgbIm->SetRegions(region);
    rgbIm->SetOrigin(origin);
    rgbIm->SetSpacing(spacing);

    vectorIm->Allocate();
    rgbIm->Allocate();

    /*Create the Iterators */
    typedef itk::ImageRegionIterator< VectorImageType >
VectorIteratorType;
    typedef itk::ImageRegionIterator< RGBImageType >    RGBIteratorType;

    VectorIteratorType    vIt(vectorIm, vectorIm->GetLargestPossibleRegion()
);
    RGBIteratorType        cIt(rgbIm, vectorIm->GetLargestPossibleRegion()
);

    /* Let's fill the vector image with random numbers between -1 and 1 */
    srand( time(NULL) );

    for( vIt.GoToBegin(); !vIt.IsAtEnd(); ++vIt ){
        VectorPixelType    setPixel;
        setPixel[0] = (VectorComponentType) rand()/RAND_MAX;  // Random
value
        setPixel[1] = (VectorComponentType) rand()/RAND_MAX;
        setPixel[2] = (VectorComponentType) rand()/RAND_MAX;

        setPixel[0] = rand() > RAND_MAX/2 ? setPixel[0] : -1 * setPixel[0];
//Randomly +'ve or -'ve
        setPixel[1] = rand() > RAND_MAX/2 ? setPixel[1] : -1 * setPixel[1];
        setPixel[2] = rand() > RAND_MAX/2 ? setPixel[2] : -1 * setPixel[2];

        vIt.Set(setPixel);
    }

    /* Here you might need to find the maximum and minimum pixel values
     * and then normalize the vector pixel components by these values.
     * Presumably your vectors are unit vectors so it is likely that the max
and
     * min values are -1 and 1, but I don't know. */

    for( cIt.GoToBegin(), vIt.GoToBegin(); !vIt.IsAtEnd(); ++cIt, ++vIt ){
        VectorPixelType getPixel;
        getPixel = vIt.Get();

        RGBPixelType rgbPixel; // since abs(values) of vector are 0-1,
multiply by RGB pixel max
        rgbPixel[0] = std::numeric_limits< RGBComponentType >::max() * fabs
(getPixel[0]);
        rgbPixel[1] = std::numeric_limits< RGBComponentType >::max() * fabs
(getPixel[1]);
        rgbPixel[2] = std::numeric_limits< RGBComponentType >::max() * fabs
(getPixel[2]);

        cIt.Set(rgbPixel);
    }

    typedef itk::ImageFileWriter< RGBImageType >    RGBWriterType;
    RGBWriterType::Pointer    writer = RGBWriterType::New();
    writer->SetInput( rgbIm );
    std::string fileName = argv[1];
    writer->SetFileName( fileName );

    writer->Update();

    return 0;
}


On Mon, Feb 14, 2011 at 8:42 AM, alba garin <albagarin1986 at hotmail.com>wrote:

>
> Hello Seth
>
> the error is "Unhandled exception at 0x012cb166 in
> tesImageReaderExample.exe: 0xC0000094: Integer division by zero." now in
> doing like this
> <code>
> rgb.SetRed(x);
> rgb.SetGreen(y);
>  rgb.SetBlue(z);
> indexIt=itRGB.GetIndex();
> rgbImage->SetPixel(indexIt,rgb);
> </code>
> but first i tried like this
> <code>
> rgb.SetRed(x); rgb.SetGreen(y); rgb.SetBlue(z);
> itRGB.Set(rgb);
> </code>
> where rgb is itk::RGBPixel<double>, x,y and z are the values extracted from
> the eigenvector, itRGB is the iterator that goes through the
> itk::Image<RGBPixelType,3>.
> thanks,
>
>
>
>
> ------------------------------
> From: seth at mech.ubc.ca
> Date: Mon, 14 Feb 2011 08:14:57 -0800
>
> Subject: Re: [Insight-users] map the colours from the eigenvector
> To: albagarin1986 at hotmail.com
> CC: insight-users at itk.org
>
>
> Hi alba,
> What's the error?
>
> Seth
>
> On Mon, Feb 14, 2011 at 3:09 AM, alba garin <albagarin1986 at hotmail.com>wrote:
>
>  Hi again,
>
> the problem is that i got an error at doing itRGB.Set(rgbPixel).
> why is that happening?
> Thanks,
>
>
>
>
>
> ------------------------------
> From: albagarin1986 at hotmail.com
> To: seth at mech.ubc.ca; insight-users at itk.org
> Date: Mon, 14 Feb 2011 07:52:29 +0000
> Subject: Re: [Insight-users] map the colours from the eigenvector
>
>
>
> Thanks, I needed to save each component in RGB pixels so your pseudo code
> is useful.
>
>
>
>
>
> ------------------------------
> From: seth at mech.ubc.ca
> Date: Fri, 11 Feb 2011 13:25:25 -0800
> Subject: [Insight-users] map the colours from the eigenvector
> To: insight-users at itk.org; albagarin1986 at hotmail.com
>
> Hi Alba,
> I'm not sure if I understand, is it the magnitude that you want mapped or
> the direction?  You can use ParaView to create visualizations of both from
> an image of vectors, and that my be the best way to go.  ParaView makes some
> very pretty pictures!
>
> If you want to save each component of the vector to the components of an
> RGB pixel, you could use iterators to move through the image and save each
> component of the eigenvector to an RGB pixel as in the following pseudo code
> (more details in the Software Guide, page 738 on).
>
> RBGIteratorType itRGB(rbgImage, vectorImage->GetLargestPossibleRegion());
> ConstVectorIteratorType constItVec(vectorImage,
> vectorImage->GetLargestPossibleRegion());
>
> for ([the whole region])
> {
>    VectorImageType::PixelType currentVectorPixel = constItVec.Get();
>
>    RGBImageType::PixelType currentRGBPixel;
>    currentRGBPixel->SetBlue(currentVectorPixel[0]); // make sure your
> vector component type and pixel type match
>    currentRGBPixel->SetGreen(currentVectorPixel[1]);  // you will also want
> to normalize your vector components by the largest in the image.
>    currentRGBPixel->SetRed(currentVectorPixel[2]);
>
>    itRGB.Set(currentRGBPixel);
> }
>
>
>
> _____________________________________ Powered by www.kitware.com Visit
> other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.html Kitware offers ITK
> Training Courses, for more information visit:
> http://www.kitware.com/products/protraining.html Please keep messages
> on-topic and check the ITK FAQ at: http://www.itk.org/Wiki/ITK_FAQ Follow
> this link to subscribe/unsubscribe:
> http://www.itk.org/mailman/listinfo/insight-users
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20110214/d0232c36/attachment.htm>