Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

I'm in the process of learning on how to remove clouds and cloud-shadows from Landsat images using GRASS. So far I was able to remove a lot of clouds from those images using i.landsat.acca, but it seems like they still require a lot of refinements:

Removing Landsat 7 clouds

I followed the GRASS manual for i.landsat.acca to obtain the image above: http://grass.osgeo.org/manuals/html70_user/i.landsat.acca.html.

Basically, the process is to run different Landsat band images for the same scene with i.landsat.toar, and then processed with i.landsat.acca to get the cloud areas. Then use the cloud area as a mask using r.mapcalc:

r.mapcalc "MASK = if(isnull(clouds))"

Problem#1: the cloud edges is still visible. I need to increase the cloud-mask buffer to catch those renegade cloud edges.

Problem#2: the cloud-shadows is still visible. If I include the -s flag to include a category for cloud-shadow when running i.landsat.acca, a large amounts of shadow from the terrain will consequently be removed as well... So to catch the cloud-shadows, I believe I need to make a new mask for them - by offsetting the current cloud mask until they cover the cloud-shadow areas.

To achieve both objectives, I believe I need to modify and play with the r.mapcalc command, but I can't figure out how to write it (due to my bad Math and my-even-worse programming skills)..

So, my ultimate question:

Question#1: How do I edit r.map.calc to increase the mask area (let's say by a buffer of 10 units)? Question#2: How do we offset a mask in r.map.calc (in this case a little to the south-west direction)?

If I know the answer to both of these: I can combine both cloud-mask and cloud-shadow-mask, to get a relatively cloud-free image. Then it's a matter of overlaying it to images of the same scene of different dates.

I appreciate any help that I can get. Many thanks!

share|improve this question
3  
Did you look at grass.osgeo.org/manuals/html70_user/r.grow.distance.html? –  markusN Oct 22 '12 at 20:27
2  
Note that the amount of offset must depend on the altitudes of the clouds and the elevation of the terrain: it is not a simple lateral movement of the clouds to mask the shadows! –  whuber Oct 22 '12 at 20:29
    
@markusN I'm not quite clear on the use of r.grow.distance.. Can we use it to offset stuff around? What should I use as the distance parameter? –  Haziq Oct 22 '12 at 21:02
1  
@whuber I've looked around for methods to do this - one of them is fs.fed.us/global/iitf/pubs/iitf-gtr32.pdf.. They've stated that "the effects of Sun illumination angle, topography, and cloud elevation should be considered"; but they didn't exactly tell how they use these parameters.. Do you know any good documentation that explains this kind of operation? –  Haziq Oct 22 '12 at 21:14
    
@whuber Oh.. They've picked shadows located over low-elevation landscapes where clouds are at the greatest height above land, project the cloud mask over it and applied a 10px buffer to account for possible variations in the shadowed areas related to variation in cloud elevation and topography.. Now I need to learn how to actually move masks around using GRASS! –  Haziq Oct 22 '12 at 21:20

1 Answer 1

up vote 3 down vote accepted

Thanks to @markusN and @whuber for pointing out some of my glaring oversight on my approach of dealing with this problem.

At first I thought what I needed was a mathematical r.mapcalc formula that could magically offset a mask around - which is the wrong thought-process. What I really need is to first apply specific raster commands to the cloud raster, and then include the processed raster in the r.mapcalc as a mask.

So, the answer to my problems:

Problem#1: the cloud edges is still visible. I need to increase the cloud-mask buffer to catch those renegade cloud edges.

Solution: Use the cloud-cover generated from i.landsat.acca and apply r.buffer or r.grow to it. For example:

r.buffer input=clouds output=clouds_buffered distances=30 unit=meters

Problem#2: the cloud-shadows is still visible. If I include the -s flag to include a category for cloud-shadow when running i.landsat.acca, a large amounts of shadow from the terrain will consequently be removed as well... So to catch the cloud-shadows, I believe I need to make a new mask for them - by offsetting the current cloud mask until they cover the cloud-shadow areas.

Solution: Duplicate the cloud-cover and offset it to its shadow using r.region:

r.region map=shadow n=n+100 e=e+100 w=w+100 s=s+100

or just by using 2-directions with + and - sign:

r.region map=shadow n=n-100 e=e+100

The final mask:

We can finally combine both cloud and cloud-shadow cover using r.patch:

r.patch input=clouds,shadow output=totalcloudcover

and then create a mask from that with r.mapcalc:

r.mapcalc "MASK = if(isnull(totalcloudcover))"

Final thoughts:

As @whuber pointed out, it's not as simple as reprojecting clouds around to mask its own shadow, because of the range of terrain, cloud-elevation and illumination angle effects.. The mask resulting from lateral movements overestimates all the areas with potential cloud shadows, and could include some underlying usable data - and should be further differentiated. At the very least, this technique is a start.

share|improve this answer
1  
+1 This reply is exceptionally thoughtful and well reasoned. Wonderful! –  whuber Oct 23 '12 at 13:48

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.