Location>code7788 >text

IDL de-cloud processing based on Landsat QA bands [code

Popularity：768 ℃/2024-08-01 14:29:28

IDL de-cloud processing based on Landsat QA bands [code

The Landsat QA band (Quality Assessment Band) is a special band in Landsat satellite image data that is present in every Landsat 5-9 product. Although we commonly use Landsat image data in B1-B7 bands, QA band is not one of them. It reflects pixels in the categories of clouds, cloud shadows, snow, etc. and is often applied in image processing for cloud pixel removal.

Recently, I have been writing about landsat pixel de-cloud processing, and after checking many QA band value explanation notes on the Internet, I found that they are all based on thebinary system (math)But IDL is different from GEE's algorithm, there is no >> this kind of bitwise operator, can only be converted to binary, and then write their own algorithm to deal with. Algorithm written, in order to send the blog went to check the official website, and found that the official website has been updated to explain the QA wave value has been updated to thedecimal system, so wrote a bit more about de-cloud processing based on decimalization (for crying out loud really).

Method 1: Processing based on the interpretation of binary values given by QA

The image above lists the meaning represented by each bit of the QA band, which is binary stored information.

The QA bands are stored in decimal, so they are converted to binary values for judgment, as illustrated below for a pixel binary value. There is a high probability that the pixel is a cloud.

Code Thoughts:

  1. Reads an image and converts the data from decimal to binary format
  2. Cloud pixels are identified and labeled, e.g. (only cloud pixels and cloud shadow pixels are removed), for convenience, only the two with bits 3 and 4 are used as reference, and no confidence values (CONFIDENCE) have been added.
  3. Create a mask to mask the original image
PRO LANDSAT_MASK_CLOUD
  COMPILE_OPT IDL2
  e = ENVI()
  
  raster = ('F:\gbytemp\LC09_L2SP_127031_20220509_20220511_02_T1\LC09_L2SP_127031_20220509_20220511_02_T1_MTL.xml')
  qaPixelRaster = ('F:\gbytemp\LC09_L2SP_127031_20220509_20220511_02_T1\LC09_L2SP_127031_20220509_20220511_02_T1_QA_PIXEL.TIF')
  data = ()
  dimensions = SIZE(data, /DIMENSIONS)
  dataBit = ()
; QA Bit Description values
; 0 Fill
; 1 Dilated Cloud 1
; 2 Cirrus 1
; 3 Cloud 1
; 4 Cloud Shadow 1
; 5 Snow 1
; 8-9 Cloud Confidence 01Low 10Reserved 11 High
; 10-11 Cloud Shadow Confidence 01Low 10Reserved 11 High
; 12-14 Snow/Ice Confidence 01Low 10Reserved 11 High
; 14-15 Cirrus Confidence 01Low 10Reserved 11 High

  stop
  mask = MAKE_ARRAY(dimensions, VALUE=1, /INTEGER)
  FOR N = 0, dimensions[0]-1 DO BEGIN
    FOR M = 0, dimensions[1]-1 DO BEGIN
; This article only uses thebit 3（surname Yun）、bit 5（surname Yun阴影）进行去surname Yun操作
; included among these3cap (a poem)4Indicates the binary position，Counting from right to left（0commencement）the reason why3cap (a poem)4The index position of the-4cap (a poem)-5
      IF dataBit[-4, N, M] EQ 1 OR dataBit[-5, N, M] EQ 1 THEN BEGIN
        mask[N, M] = 0
      ENDIF
    ENDFOR
  ENDFOR
  

  file = ()
  maskRaster = ENVIRaster(mask, URI=file)
  
  maskedRaster = ENVIMaskRaster(raster[0], maskRaster)

  , maskedRaster
  view=()
  layer=(maskedRaster)
  stop
END

​ Comparison chart of de-clouded results:

Method 2: Processing based on the interpretation of decimal values given by QA

The decimal values explain the meaning as follows:

Code Thoughts:

  1. Read image
  2. Cloud pixels are identified, and labeled, e.g. (only cloud pixels and cloud shadow pixels are removed), and for convenience, only the high confidence cloud 22280, and the high confidence cloud shadow 23888 are used as references, the
  3. Create a mask to mask the original image
PRO LANDSAT_MASK_CLOUD
  COMPILE_OPT IDL2
  e = ENVI()
  
  raster = ('F:\gbytemp\LC09_L2SP_127031_20220509_20220511_02_T1\LC09_L2SP_127031_20220509_20220511_02_T1_MTL.xml')
  qaPixelRaster = ('F:\gbytemp\LC09_L2SP_127031_20220509_20220511_02_T1\LC09_L2SP_127031_20220509_20220511_02_T1_QA_PIXEL.TIF')
  data = ()
  dimensions = SIZE(data, /DIMENSIONS)

  stop
  mask = MAKE_ARRAY(dimensions, VALUE=1, /INTEGER)
  FOR N = 0, dimensions[0]-1 DO BEGIN
    FOR M = 0, dimensions[1]-1 DO BEGIN
      IF data[N, M] EQ 55052 OR data[N, M] EQ 23888 THEN BEGIN
        mask[N, M] = 0
      ENDIF
    ENDFOR
  ENDFOR

  file = ()
  maskRaster = ENVIRaster(mask, URI=file)
  
  maskedRaster = ENVIMaskRaster(raster[0], maskRaster)

  , maskedRaster
  view=()
  layer=(maskedRaster)
  stop
END

​ Comparison chart of de-clouded results:

<Previous Next>