How to create the matrix for every image pixel from HSI image in matlab?

Hello guys,

I have Hyperspectral imaging file (X x Y x Z) ,How to create the matrix for every image pixel from HSI image in matlab?

The HSI file size is bigger than 5 MB (Matlab attachment limit), so i have uploaded the file on drive Link to access image file: https://drive.google.com/drive/folders/1GTSkQMog4EvffH_geeKJun0oX6BCh5zo?usp=sharing

any idea thanks a ton

Emma

[Emma Stewart]

@Image Analyst, do you have any idea , thanks a lot for help. i saw your answers and i think, you may knows the way to do it.

Answer 1

[Walter Roberson]

https://www.mathworks.com/matlabcentral/answers/853385-reading-a-raw-file-from-a-hyperspectral-camera-using-envi-information#answer_722390

[Emma Stewart]

@Walter Roberson, thanks a lot. I read the answer but i would like to create matrix for every pixel and that code gives the only information height X width X number of bands. Can you please tell me how to do it further?

Sorry, if i did not explained it very well

Emma

[Image Analyst]

How many elements would each pixel's "matrix" hold? What would be the dimensions of that matrix?

[Emma Stewart]

@Image Analyst , the image has dimensions (353 320 256) =(X*Y*Number of bands) , so it should carry all pixels. I just would like to use this matrix again for clustering purpose.

Thanks a lot

Emma

[Image Analyst]

OK, so those are the image dimensions. And you have

numPixels = 353 * 320 * 256

numPixels = 28917760

pixels. And you want a matrix for every of those 28.9 million pixels. What would one of those 28.9 million matrices look like? Each pixel has a matrix (you said) so what is in the matrix? And how do you want it organized? The matrices would also have to be in a matrix because you can't have 28.9 million differently named variables.

And a hyperspectral image has 4 or more wavelength channel images. So you can't get an RGB or HSI or HSV image from that unless you do some kind of combining or selection of the channels to get down to only 3 channels.

[Emma Stewart]

@Image Analyst I would like to use this matrix for classfication purpose. (The dimension itself (X=353, Y=320, Number of bands=256). So the pixel are (X=353, Y=320) and every pixel should represent the number of bands. So i can use it that matrix for classfication purpose.

is it understandable ? or can i explain it another way?

Thanks a lot

Emma

[Image Analyst]

OK, so you don't want a matrix for every pixel. You want a classified image that gives the class of each pixel in the image. Well how many classes to you think you have? Have you tried kmeans or a discriminant classifier or KNN? I'm attaching demos for each.

[Emma Stewart]

@Image Analyst , I want matrix for each pixel. I just want to use that matrix for clustering purpose. Can you please tell me how to do it in the my attached image

[Walter Roberson]

num2cell passing in dimension number of 3 as the second parameter will give you a cell array with one entry per pixel.

That is what you are asking for. It will also be useless for you, as you cannot train one pixel at a time.

I would suggest that what you more likely want is to reshape(data, [], size(data, 3)) (which might possibly need to be transposed afterwards)

[Emma Stewart]

@Walter Roberson I did not understand. I meant not one pixel but all the pixel values of image. Then use that pixel band values for clustering purpose. (See image)

[Walter Roberson]

If you wanted to split into bands, equivalent to

for K = 1 : size(data,3)
    parts{K} = data(:,:,K);
end

then you can just use the loop, or you can use num2cell() with an appropriate second argument.

[Walter Roberson]

data = rand(3,4,5);
p = num2cell(data,1); size(p), size(p{1})

ans = 1×3

1 4 5

ans = 1×2

3 1

p = num2cell(data,2); size(p), size(p{1})

ans = 1×3

3 1 5

ans = 1×2

1 4

p = num2cell(data,3); size(p), size(p{1})

ans = 1×2

3 4

ans = 1×3

1 1 5

p = num2cell(data, [1 2]); size(p), size(p{1})

ans = 1×3

1 1 5

ans = 1×2

3 4

[Emma Stewart]

@Walter Roberson , I did not understand, p is the pixel output ?

[Walter Roberson]

The above illustrates that you can use num2cell to split up your array into cell arrays along any subset of the dimensions.

num2cell(ARRAY, DIMENSIONLIST) is going to split up ARRAY into cells, giving a cell array. The size() of the cells returned will be size(ARRAY) with each dimension mentioned in DIMENSIONLIST replaced by 1. The size() of each entry inside the cell array will be size(ARRAY) with each dimension listed in DIMENSIONLIST preserved and the other dimensions all replaced by 1.

Hence, num2cell(data, [1 2]) results in a 1 x 1 x size(data,3) cell array, and the cells are size(data,1) by size(data,2) where those 1 and 2 are the numbers from the [1 2].

Use this to split your array into cells however is appropriate for you.

Note: you might want to cellfun(@squeeze, p) afterwards, to get rid of the singular dimensions.

[Emma Stewart]

@Walter Roberson how can i use this data for dimensionality reduction purpose? I tried but mind stuck up to imagine. Thank you in advance.

Emma

[Walter Roberson]

@Emma Stewart

I do not know. I am lost about what you are trying to do. The things that I thought you wanted to do, you would normally proceed by reshape(data, [], size(data, 3)) like I suggested yesterday, but which you did not feel was appropriate.

The general approach is that you would figure out which parts of the data should be considered to be independent samples, and which parts are considered to be the readings for each sample. You would turn that into a 2d array, in which each row was an independent sample, and each column was a reading (variable) from that sample. Then you would use pca() pca

My understanding is that you want to consider the set of spectra for each pixel as the readings, and you want to consider each pixel as independent, with order of the pixels within the image not relevant. Such a situation might occur, for example, if you had hyperspectral satellite readings and wanted to classify each pixel as (say) ground, or normal water, or acid water: the position of the pixel within the image might help hint at classification but is not important for training purposes.

If that was indeed the kind of task you are trying to do, where position is not part of the training information but the spectra is, then reshape(data, [], size(data, 3))

[Emma Stewart]

@Walter Roberson , Thanks a lot for explaining. Yes i want to classify each pixel as different class (ground, or normal water, or acid water). I just want to reduce the dimesionality reduction , so i can use this data of pixel for classification (The position of the pixel within the image might help hint at classification)

Thanks a lot

Emma

[Walter Roberson]

Then you will want to do the reshape() I indicated earlier.

https://www.mathworks.com/help/stats/dimensionality-reduction.html

[Emma Stewart]

@Walter Roberson can you please write code here. I tried but unable to obtain the cluster image. I am totally confuse now

Thanks a lot.

Emma

[Walter Roberson]

Simple approach

numclasses = 3;
M = reshape(Your3Ddata, [], size(data, 3));
[idx, c] = kmeans(M, num_classes+1);

xcol = 2;
ycol = 7;
gmatrixplot(M(:,xcol), M(:, ycol), idx);

xcol and ycol are selected data columns for plotting purposes.

This example does not attempt any dimension reduction.

You might notice that I clustered with one more cluster than the number of data groups. It is recommended to have one extra cluster for "other".

Note that kmeans by itself is not able to tell you what the clusters "mean" and the cluster centers will not be in any particular order.

You might want to normalize the data per column before clustering, unless you want one column to be more important than another (the columns you want to be more important should have larger variance)

[Esila Darci]

@Walter Roberson@Image Analyst @Emma Stewart

I was following this interesting insightful discussion and i have my doubt with respect to image attached in the question. How about the , if I just want this matrix (Number of Pixels X wavelength value), which is (112960 x 256)

Here number of pixels is (353 X 320 = 112960) and

respective wavelength value at every 256 bands, which is in between (917.32 to 1742.34).

How to create a code for to generate above matrix?

Thanks a lot for help.

[Esila Darci]

@Walter Roberson@Image Analyst @Emma Stewart

1) I was following this interesting insightful discussion and i have my doubt with respect to image attached in the question. How about the , if I just want this matrix (Number of Pixels X wavelength value), which is (112960 x 256)

Here number of pixels is (353 X 320 = 112960) and

respective wavelength value at every 256 bands, which is in between (917.32 to 1742.34).

How to create a code for to generate above matrix?

Thanks a lot for help.

2) Also for me this file is not readable with code you mentioned earlier email @Walter Roberson? (link you mentioned: Reading a RAW file from a hyperspectral camera using ENVI information - (mathworks.com) I do not know why, this error shows to while open the image mentioned in the question.

Error using enviinfo>getKeyValuePair

Invalid expression in line 5 of input file.

Error in enviinfo (line 133)

[key, value] = getKeyValuePair(lines);

Error in wavelengthperpixel (line 2)

info = enviinfo('Labsample.hdr');

Any idea in change in code would be highly apreciable, thank you

[Walter Roberson]

LabsampleImage = ....
data = reshape(LabsampleImage, [], size(LabsampleImage,3));

Now process data

[Walter Roberson]

If I recall correctly, there has been more than one version of ENVI standards; it is possible that the one you are using is not supported by enviinfo()

Can you attach the .hdr file? zip it and attach the zip .

[Esila Darci]

@Walter Roberson I have attached .hdr file link here. The file size is more than 5 MB so i have to upload it on google drive

HSI image - Google Drive

Can you please try to open this one.

Thanks

[Image Analyst]

@Esila Darci:

Can you start your own thread, rather than hijack @Emma Stewart's thread? And say how many classes you want to find in your hyuper spectral image. It can get confusing when giving answers to two different people's problems in the same discussion.

[Esila Darci]

@Image Analyst Apolosige for consfusing you.

yes , i can start it but i thought my problem is small part of the question asked by Emma.

Also my question is not about the class but this one

1) I have my own doubt with respect to image attached in the question. How about the , if I just want this matrix (Number of Pixels X wavelength value), which is (112960 x 256)

Here number of pixels is (353 X 320 = 112960) and

respective wavelength value at every 256 bands, which is in between (917.32 to 1742.34).

How to create a code for to generate above matrix?

Thanks a lot for help.

2) Also for me this file is not readable with code you mentioned earlier comment by walter robbinson. (link mentioned: Reading a RAW file from a hyperspectral camera using ENVI information - (mathworks.com) I do not know why, this error shows to while open the image mentioned in the question.

Error using enviinfo>getKeyValuePair

Invalid expression in line 5 of input file.

Error in enviinfo (line 133)

[key, value] = getKeyValuePair(lines);

Error in wavelengthperpixel (line 2)

info = enviinfo('Labsample.hdr');

Many thanks in advance

[Walter Roberson]

@Esila Darci

If you look on line 5 of your .hdr file you will find a "sensor type =" line that gives a software version. But then if you look on the next line, there is a DAQ version with no keyword and "=" at the beginning. The version information has been split across two lines, and that is a problem. You need to edit the .hdr file and join those two lines into a single line.

[Esila Darci]

@Walter Roberson, i tried but error is the same. I have did this , am i right? can you please attach your screenshot as well, so i can do it same like you.

Thanks a lot

[Walter Roberson]

ENVI
description = {
File Imported into ENVI }
file type = ENVI
sensor type = Specim Hyperspectral Sensor, SW ver 3.6219525 DAQ SW ver 3.62.195.25-R
acquisition date = DATE(dd-mm-yyyy): 07-04-2022
qpfTiming = {  13.43, 8273.93}

and then I used

info = enviinfo('Sample.hdr');
data = multibandread('Sample.raw',[info.Height, info.Width, info.Bands],...
           info.DataType, info.HeaderOffset, info.Interleave, info.ByteOrder);

[Esila Darci]

@Walter Roberson

Thanks a lot.

Now i got it,

What about the Ieee-ie or Ieee-be should i use in this case ? I know the sequence of byte will change but what should i prefer in my case, do you have any idea.

Also what should i select in between this three option in my case

• 'bsq' — Band-Sequential
• 'bil'— Band-Interleaved-by-Line
• 'bip'— Band-Interleaved-by-Pixel

Again, thanks a lot for help..

Esila

[Walter Roberson]

The .hdr file itself contains the information about the interleave; it is 'bil' for that file. enviinfo() retrieves that information and stores it in the info structure, and the next line of the code passes the information info.Interleave to the multibandread() call.

Likewise, the byte order information is in the header, which says 0; that is the info.ByteOrder that will be passed to the call.

[Esila Darci]

@Walter Roberson many thanks for detailed response. One quick question as you mentioned " Likewise, the byte order information is in the header, which says 0; that is the info.ByteOrder that will be passed to the call."

So it is "Ieee-ie" or "Ieee-be" should i use, if i would like to extract information with manual function?

Thanks a lot

Esila

[Walter Roberson]

https://www.l3harrisgeospatial.com/docs/enviheaderfiles.html

0 corresponds to least significant byte first, which is ieee-le

[Esila Darci]

Thanks @Walter Roberson