Help with code for hyperspectral image analysis

I have BIL image data of an estuary with 1674 lines, 1302 samples, and 156 bands. I need to mask everything that is land, which is easily identified by its high reflectance values in certain bands - I use band 140 to identify it in my code.

Here's the code that I've worked up, and it does indeed write an ENVI file, but it's reduced to a single band instead of having all 156 in the unmasked areas. How do I get my output file to have all 156 bands??

Parts of this code were written by someone else, so just because I did something once in this code that doesn't mean that I understand it or could apply the concept again!

Thank you so much!
-K


pro envi_masking_image

;select your image file

;Get your file
image=dialog_pickfile()

;check to see if valid for ENVI
envi_open_file, image, r_fid=fid

;learn about the image
ENVI_FILE_QUERY, fid, ns=ns, nl=nl, nb=nb, data_type=data_type
help, data_type

;build array based on data type
data_holder=intarr(ns, nb, nl)
help, data_holder

;open a new file for reading as IDL file unit number 1
openr, 1, image

;write the data in my image file to the file 'data_array'
readu, 1, data_holder
close, 1

;create the mask
;locate the place that you want data from.
;-1 is just because, then the next two numbers specify the first sample, last sample, first line, and last line
;that you will read from. Envi starts with 0, so to read the whole image read from 0 to ns-1 and nl-1.
dims = [-1, 0, ns-1, 0, nl-1]

;get data from the image, pos seems to indicate which band to get data from
f = ENVI_GET_DATA(fid=fid, dims=dims, pos=140)

lcmask = f LT 40


;apply the mask to the image array

masked_image=convol(data_holder, lcmask)

;write out the result as an envi file
envi_write_envi_file, masked_image, out_name='masktest'


end

Your problem is your call to ENVI_GET_DATA. POS has to be a scalar, so ENVI_GET_DATA returns only a single band at a time. You should loop through all the bands, apply your filter, then write the band to the output file using writeu. Then, call ENVI_SETUP_HEAD to write the header file.

-B

On Monday, June 11, 2012 10:59:55 AM UTC-4, Kel wrote:

> ;get data from the image, pos seems to indicate which band to get data from
> f = ENVI_GET_DATA(fid=fid, dims=dims, pos=140)
>
> lcmask = f LT 40
>
>
> ;apply the mask to the image array
>
> masked_image=convol(data_holder, lcmask)
>
> ;write out the result as an envi file
> envi_write_envi_file, masked_image, out_name='masktest'
>
>
> end

Thanks for your reply! I'm sorry to say that I still don't understand, though.

Calling up a single band with ENVI_GET_DATA was intentional, because I wanted to make a single-band file showing ones for locations where that band was less than 40, and zeros for where it's greater than 40.

The part that I'm stuck on is how I apply that information to all pixels inmy image. At each pixel, if the band value at 140 is greater than 40, I want all 156 bands zeroed out. If it's less than 40 at that one band, then itshould be left alone.

Any further explanation you can give would be very appreciated!



On Monday, June 11, 2012 12:43:00 PM UTC-4, Brian J. Daniel wrote:
> Your problem is your call to ENVI_GET_DATA. POS has to be a scalar, so ENVI_GET_DATA returns only a single band at a time. You should loop throughall the bands, apply your filter, then write the band to the output file using writeu. Then, call ENVI_SETUP_HEAD to write the header file.
>
> -B
>
> On Monday, June 11, 2012 10:59:55 AM UTC-4, Kel wrote:
>
> > ;get data from the image, pos seems to indicate which band to get data from
> > f = ENVI_GET_DATA(fid=fid, dims=dims, pos=140)
> >
> > lcmask = f LT 40
> >
> >
> > ;apply the mask to the image array
> >
> > masked_image=convol(data_holder, lcmask)
> >
> > ;write out the result as an envi file
> > envi_write_envi_file, masked_image, out_name='masktest'
> >
> >
> > end

Try this, I didn't test it...
Cheers, Klemen

pro envi_masking_image

;select your image file

;Get your file
image=dialog_pickfile()

;check to see if valid for ENVI
envi_open_file, image, r_fid=fid

;learn about the image
ENVI_FILE_QUERY, fid, ns=ns, nl=nl, nb=nb, data_type=data_type
help, data_type

;build array based on data type
data_holder=intarr(ns, nb, nl)

;make a mask
f = ENVI_GET_DATA(fid=fid, dims=dims, pos=140)
indx = where(f LT 40, count)

;applay the mask to all bands
if count gt 0 then begin
for i=0,nb-1 do begin
f = ENVI_GET_DATA(fid=fid, dims=dims, pos=i
f[indx] = 0
data_holder[*,i,*] = f
endfor
endif

;write out the result as an envi file
envi_write_envi_file, data_holder, out_name='masktest'


end

See my adjustments to Klemen's code below. This way you don't have to load the entire image, masked or otherwise, into memory. Also, untested.

pro envi_masking_image

;select your image file

;Get your file
image=dialog_pickfile()

; Pick an output filename
output = dialog_pickfile()

;check to see if valid for ENVI
envi_open_file, image, r_fid=fid
IF fid[0] EQ -1 THEN Message,'Invalid file name: '+image

;learn about the image
ENVI_FILE_QUERY, fid, ns=ns, nl=nl, nb=nb, data_type=data_type,Interleave=interleave
help, data_type

;build array based on data type
;data_holder=Make_Array(ns, nb, nl, Type=data_type) ; memory intensive

;make a mask
f = ENVI_GET_DATA(fid=fid, dims=dims, pos=140)
indx = where(f LT 40, count)
f[indx] = 0 ; now f is your mask

;apply the mask to all bands and save
OpenW, lun, output, /Get_Lun
if count gt 0 then begin
for i=0,nl-1 do begin
data = ENVI_GET_DATA(fid=fid, dims=dims, pos=i)
writeu,lun,data*f

endfor
endif

; close binary file
close,lun
free_lun, lun

; write out an envi header file
envi_setup_head, Data_Type=data_type, /Open, R_FID=masked_fid, $
/WRITE, FNAME=output, Interleave=2 ; BIP

end

Thanks you both so much! I had just gotten an out-of-memory error on Klemen's code, which seemed to be going ok up until that point, so I will try this modification now and hope it helps.

On Wednesday, June 20, 2012 3:55:08 PM UTC-4, Brian J. Daniel wrote:
> See my adjustments to Klemen's code below. This way you don't have to load the entire image, masked or otherwise, into memory. Also, untested.
>
> pro envi_masking_image
>
> ;select your image file
>
> ;Get your file
> image=dialog_pickfile()
>
> ; Pick an output filename
> output = dialog_pickfile()
>
> ;check to see if valid for ENVI
> envi_open_file, image, r_fid=fid
> IF fid[0] EQ -1 THEN Message,'Invalid file name: '+image
>
> ;learn about the image
> ENVI_FILE_QUERY, fid, ns=ns, nl=nl, nb=nb, data_type=data_type,Interleave=interleave
> help, data_type
>
> ;build array based on data type
> ;data_holder=Make_Array(ns, nb, nl, Type=data_type) ; memory intensive
>
> ;make a mask
> f = ENVI_GET_DATA(fid=fid, dims=dims, pos=140)
> indx = where(f LT 40, count)
> f[indx] = 0 ; now f is your mask
>
> ;apply the mask to all bands and save
> OpenW, lun, output, /Get_Lun
> if count gt 0 then begin
> for i=0,nl-1 do begin
> data = ENVI_GET_DATA(fid=fid, dims=dims, pos=i)
> writeu,lun,data*f
>
> endfor
> endif
>
> ; close binary file
> close,lun
> free_lun, lun
>
> ; write out an envi header file
> envi_setup_head, Data_Type=data_type, /Open, R_FID=masked_fid, $
> /WRITE, FNAME=output, Interleave=2 ; BIP
>
> end

Update: Still running out of memory. I've already changed the values in the idlde.ini file, as recommended in other posts, but now I'm stuck. Google is suggesting that I increase the virtual memory of my pc, so I will give that a try.

On Wednesday, June 20, 2012 3:59:36 PM UTC-4, Kel wrote:
> Thanks you both so much! I had just gotten an out-of-memory error on Klemen's code, which seemed to be going ok up until that point, so I will try this modification now and hope it helps.
>
> On Wednesday, June 20, 2012 3:55:08 PM UTC-4, Brian J. Daniel wrote:
> > See my adjustments to Klemen's code below. This way you don't have to load the entire image, masked or otherwise, into memory. Also, untested.
> >
> > pro envi_masking_image
> >
> > ;select your image file
> >
> > ;Get your file
> > image=dialog_pickfile()
> >
> > ; Pick an output filename
> > output = dialog_pickfile()
> >
> > ;check to see if valid for ENVI
> > envi_open_file, image, r_fid=fid
> > IF fid[0] EQ -1 THEN Message,'Invalid file name: '+image
> >
> > ;learn about the image
> > ENVI_FILE_QUERY, fid, ns=ns, nl=nl, nb=nb, data_type=data_type,Interleave=interleave
> > help, data_type
> >
> > ;build array based on data type
> > ;data_holder=Make_Array(ns, nb, nl, Type=data_type) ; memory intensive
> >
> > ;make a mask
> > f = ENVI_GET_DATA(fid=fid, dims=dims, pos=140)
> > indx = where(f LT 40, count)
> > f[indx] = 0 ; now f is your mask
> >
> > ;apply the mask to all bands and save
> > OpenW, lun, output, /Get_Lun
> > if count gt 0 then begin
> > for i=0,nl-1 do begin
> > data = ENVI_GET_DATA(fid=fid, dims=dims, pos=i)
> > writeu,lun,data*f
> >
> > endfor
> > endif
> >
> > ; close binary file
> > close,lun
> > free_lun, lun
> >
> > ; write out an envi header file
> > envi_setup_head, Data_Type=data_type, /Open, R_FID=masked_fid, $
> > /WRITE, FNAME=output, Interleave=2 ; BIP
> >
> > end

I'm testing code on a tiny file now, and it's getting really close! Still not quite what I need though.
I've simplified the operations that I want to do on the image data just to get the basic code set up, because something is going weird in writing the output file. Masking high-magnitude pixels was the original goal, now I'm just trying to increase the magnitude of each band by 2.

Goal: open an ENVI image file, use a for-loop to do a simple math operation(add 2), create a new ENVI output image that looks just like the old one but with a slightly higher magnitude across all bands.

Result: The output image file displays the expected samples, lines, and bands, but the band values have no relationship to the band values of the input file, as far as I can tell. How did they get so crazy?

Thanks again for the help, I would never have managed this on my own!

;select your image file
image=filepath('tinyfile', root_dir='d:', subdirectory=['idl'])

output=filepath('tiny2', root_dir='d:', subdirectory=['idl'])

;check to see if valid for ENVI
envi_open_file, image, r_fid=fid

;learn about the image
ENVI_FILE_QUERY, fid, ns=ns, nl=nl, nb=nb, data_type=data_type, interleave=interleave, map_info=map_info
help, data_type

;build array based on data type
data_holder=intarr(ns, nb, nl)
help, data_holder

;open a new file for reading as IDL file unit number 1
openr, 1, image

;write the data in my image file to the file 'data_array'
readu, 1, data_holder
close, 1

dims = [-1, 0, ns-1, 0, nl-1]

OpenW, lun, output, /Get_Lun
for i=0,nb-1 do begin
data = ENVI_GET_DATA(fid=fid, dims=dims, pos=i)
newdata = data+2
writeu, lun, newdata

endfor

; close binary file
close, lun
free_lun, lun

; write out an envi header file
envi_setup_head, Data_Type=data_type, /Open, R_FID=masked_fid, /WRITE, FNAME=output, Interleave=1, ns=ns, nl=nl, nb=nb







On Wednesday, June 20, 2012 4:24:54 PM UTC-4, Kel wrote:
> Update: Still running out of memory. I've already changed the values in the idlde.ini file, as recommended in other posts, but now I'm stuck. Googleis suggesting that I increase the virtual memory of my pc, so I will give that a try.
>
> On Wednesday, June 20, 2012 3:59:36 PM UTC-4, Kel wrote:
> > Thanks you both so much! I had just gotten an out-of-memory error on Klemen's code, which seemed to be going ok up until that point, so I will trythis modification now and hope it helps.
> >
> > On Wednesday, June 20, 2012 3:55:08 PM UTC-4, Brian J. Daniel wrote:
> > > See my adjustments to Klemen's code below. This way you don't have to load the entire image, masked or otherwise, into memory. Also, untested.
> > >
> > > pro envi_masking_image
> > >
> > > ;select your image file
> > >
> > > ;Get your file
> > > image=dialog_pickfile()
> > >
> > > ; Pick an output filename
> > > output = dialog_pickfile()
> > >
> > > ;check to see if valid for ENVI
> > > envi_open_file, image, r_fid=fid
> > > IF fid[0] EQ -1 THEN Message,'Invalid file name: '+image
> > >
> > > ;learn about the image
> > > ENVI_FILE_QUERY, fid, ns=ns, nl=nl, nb=nb, data_type=data_type,Interleave=interleave
> > > help, data_type
> > >
> > > ;build array based on data type
> > > ;data_holder=Make_Array(ns, nb, nl, Type=data_type) ; memory intensive
> > >
> > > ;make a mask
> > > f = ENVI_GET_DATA(fid=fid, dims=dims, pos=140)
> > > indx = where(f LT 40, count)
> > > f[indx] = 0 ; now f is your mask
> > >
> > > ;apply the mask to all bands and save
> > > OpenW, lun, output, /Get_Lun
> > > if count gt 0 then begin
> > > for i=0,nl-1 do begin
> > > data = ENVI_GET_DATA(fid=fid, dims=dims, pos=i)
> > > writeu,lun,data*f
> > >
> > > endfor
> > > endif
> > >
> > > ; close binary file
> > > close,lun
> > > free_lun, lun
> > >
> > > ; write out an envi header file
> > > envi_setup_head, Data_Type=data_type, /Open, R_FID=masked_fid, $
> > > /WRITE, FNAME=output, Interleave=2 ; BIP
> > >
> > > end

Hi,

I am not so fit in the ENVI programming, but as far I see you might have a problem with the file format in envi_setup_head. Writeu writes in a for loop first the 1st band, then the 2nd... Brian suggested in his code using Interleave=2 (BIP), I guess you should have 1 (BSQ), which you actually do, so we are not getting anywhere...

But try to check if this might not be a problem.

Cheers, Klemen

That's a good thought, thanks!

I don't have a good understanding of what exactly the writeu is doing when it runs in the for-loop. The input file is in BIL format (interleave=1) and I would like the output file to be, as well. Is there a way to specify in writeu that it should be organized one way or another?


On Thursday, June 21, 2012 5:44:22 PM UTC-4, Klemen wrote:
> Hi,
>
> I am not so fit in the ENVI programming, but as far I see you might have a problem with the file format in envi_setup_head. Writeu writes in a for loop first the 1st band, then the 2nd... Brian suggested in his code using Interleave=2 (BIP), I guess you should have 1 (BSQ), which you actually do, so we are not getting anywhere...
>
> But try to check if this might not be a problem.
>
> Cheers, Klemen

Writeu performs a direct transfer with no processing of any kind being doneto the data. So if you open just one band, process it and write it down, you will write the first line, the second line.... till the last band of thefile. Then the for loop moves to next band where this is repeated.

So as far as I understand all this envi formats, this means that you write down BSQ file. This means that you should (I have just checked it in the ENVI help, see below) set the Interleave to 0 in your last line.

This is how you get your results written correctly, I thnik. To convert it to BIl, take a look at:
http://geol.hu/data/online_help/Conv...IL_BIP%29.html
But I guess that you could have a memory problem again....



INTERLEAVE
Set this keyword to one of the following integer values to specify the interleave output:

0: BSQ
Band-sequential; an interleave format where each line of the data is followed immediately by the next line in the same spectral band. This format is optimal for spatial (x,y) access of any part of a single spectral band.
1: BIL
Band-interleaved-by-line; an interleave format that stores the first line of the first band, followed by the first line of the second band, followed by the first line of the third band, interleaved up to the number of bands.
2: BIP
Band-interleaved-by-pixel; an interleave format that stores the first pixelfor all bands in sequential order, followed by the second pixel for all bands, followed by the third pixel for all bands, and so forth, interleaved up to the number of pixels.

On Friday, June 22, 2012 3:44:26 AM UTC-4, Klemen wrote:
> Writeu performs a direct transfer with no processing of any kind being done to the data. So if you open just one band, process it and write it down,you will write the first line, the second line.... till the last band of the file. Then the for loop moves to next band where this is repeated.
>
> So as far as I understand all this envi formats, this means that you write down BSQ file. This means that you should (I have just checked it in the ENVI help, see below) set the Interleave to 0 in your last line.
>
> This is how you get your results written correctly, I thnik. To convert it to BIl, take a look at:
> http://geol.hu/data/online_help/Conv...IL_BIP%29.html
> But I guess that you could have a memory problem again....
>
>
>
> INTERLEAVE
> Set this keyword to one of the following integer values to specify the interleave output:
>
> 0: BSQ
> Band-sequential; an interleave format where each line of the data is followed immediately by the next line in the same spectral band. This format isoptimal for spatial (x,y) access of any part of a single spectral band.
> 1: BIL
> Band-interleaved-by-line; an interleave format that stores the first lineof the first band, followed by the first line of the second band, followedby the first line of the third band, interleaved up to the number of bands..
> 2: BIP
> Band-interleaved-by-pixel; an interleave format that stores the first pixel for all bands in sequential order, followed by the second pixel for all bands, followed by the third pixel for all bands, and so forth, interleavedup to the number of pixels.

Good point, Klemen. You should use interleave of BSQ. Once you have the ENVI file written in BSQ, use CONVERT_INPLACE_DOIT to convert the file to any interleave you like.

-Brian

That did it! Setting the interleave on the output to 0 for BSQ gave me an image with all the right spectral profiles, and I was then able to add the masking code back in to accomplish my original goal.

Thank you both so much!


On Friday, June 22, 2012 7:59:12 AM UTC-4, Brian J. Daniel wrote:
> On Friday, June 22, 2012 3:44:26 AM UTC-4, Klemen wrote:
> > Writeu performs a direct transfer with no processing of any kind being done to the data. So if you open just one band, process it and write it down, you will write the first line, the second line.... till the last band ofthe file. Then the for loop moves to next band where this is repeated.
> >
> > So as far as I understand all this envi formats, this means that you write down BSQ file. This means that you should (I have just checked it in the ENVI help, see below) set the Interleave to 0 in your last line.
> >
> > This is how you get your results written correctly, I thnik. To convertit to BIl, take a look at:
> > http://geol.hu/data/online_help/Conv...IL_BIP%29.html
> > But I guess that you could have a memory problem again....
> >
> >
> >
> > INTERLEAVE
> > Set this keyword to one of the following integer values to specify the interleave output:
> >
> > 0: BSQ
> > Band-sequential; an interleave format where each line of the data is followed immediately by the next line in the same spectral band. This format is optimal for spatial (x,y) access of any part of a single spectral band.
> > 1: BIL
> > Band-interleaved-by-line; an interleave format that stores the first line of the first band, followed by the first line of the second band, followed by the first line of the third band, interleaved up to the number of bands.
> > 2: BIP
> > Band-interleaved-by-pixel; an interleave format that stores the first pixel for all bands in sequential order, followed by the second pixel for all bands, followed by the third pixel for all bands, and so forth, interleaved up to the number of pixels.
>
> Good point, Klemen. You should use interleave of BSQ. Once you have theENVI file written in BSQ, use CONVERT_INPLACE_DOIT to convert the file to any interleave you like.
>
> -Brian