[GRASS-user] question on Stray light Landsat 8 data

Markus Metz markus.metz.giswork at gmail.com
Tue Aug 20 13:13:36 PDT 2019


Dear Gabriel,

On Tue, Aug 20, 2019 at 9:50 PM Gabriel Cotlier <gabiklm01 at gmail.com> wrote:
>
> Dear Markus,
> So I run the command g.region w=179:59:45W e=180:00:15E
>
> Now my log now before and after the command g.region w=179:59:45W
e=180:00:15E is as follows:
>
> 360 degree EW extent is exceeded by 1 cells
> 360 degree EW extent is exceeded by 1 cells
> (Tue Aug 20 16:41:45 2019) Command finished (0 sec)

> (Tue Aug 20 16:43:47 2019)

> g.region w=179:59:45W e=180:00:15E

> 360 degree EW extent is exceeded by 1.99983 cells
> 360 degree EW extent is exceeded by 1 cells
> (Tue Aug 20 16:43:48 2019) Command finished (0 sec)

> (Tue Aug 20 16:44:18 2019)

> g.region -p

> projection: 3 (Latitude-Longitude)
> zone:       0
> datum:      wgs84
> ellipsoid:  wgs84
> north:      75:00:15N
> south:      65:00:15S
> west:       179:59:45W
> east:       180:00:15E
> nsres:      0:00:30
> ewres:      0:00:30
> rows:       16801
> cols:       43200
> cells:      725803200
> 360 degree EW extent is exceeded by 1.99983 cells
> (Tue Aug 20 16:44:18 2019) Command finished (0 sec)

With GRASS 7.6, I can not reproduce this message.
g.region -p n=75:00:15N s=65:00:15S w=179:59:45W e=180:00:15E res=00:00:30

gives me
projection: 3 (Latitude-Longitude)
zone:       0
datum:      wgs84
ellipsoid:  wgs84
north:      75:00:15N
south:      65:00:15S
west:       179:59:45W
east:       180:00:15E
nsres:      0:00:30
ewres:      0:00:30
rows:       16801
cols:       43200
cells:      725803200

without any messages that 360 degree EW extent is exceeded. Which GASS
version are you using? I tested with GRASS 7.6 and GRASS 7.9.

Markus M

>
> Now appears to say that is exceeded by 1.99983 cells.... why could this
be happening?
> Thanks a lot
>
> Regards
> Gabriel
>
>
>
>
>
>
>
>
>
>
>
>
> On Tue, Aug 20, 2019 at 4:33 PM Gabriel Cotlier <gabiklm01 at gmail.com>
wrote:
>>
>> Dear Markus,
>> Thanks a lot for your response and explanation.
>> Changing the region to w=179:59:45W e=180:00:15E, am I not only avoiding
the warning, but also changing the layers to be physically correct, right?
>>
>> Thanks again for your help.
>> regards,
>> Gabriel
>>
>>
>> On Tue, Aug 20, 2019 at 4:24 PM Markus Metz <
markus.metz.giswork at gmail.com> wrote:
>>>
>>> Dear Gabriel,
>>>
>>> On Tue, Aug 20, 2019 at 12:19 AM Gabriel Cotlier <gabiklm01 at gmail.com>
wrote:
>>> >
>>> > Dear Markus,
>>> >
>>> > Thanks a lot for the clarification and explanation, your response was
indeed helpful.
>>> >
>>> > I got for all maps in the mapset I used, for both the DMSP original
raster layers and the intercallibrated rasrer layers the following:
>>> >
>>> > r.info map= name_of_raster_map
>>> >
>>> > 360 degree EW extent is exceeded by 1 cells
>>> > 360 degree EW extent is exceeded by 1 cells
>>> >
>>> > Which, following what you said before in your response I understand
makes it correct region, right?
>>>
>>> this region is correct considering the resolution with is now exactly
30 arc seconds.
>>>
>>> this region is not correct considering that 360 degree EW extent is
exceeded by 1 cell. The first column from 180:00:15W to 179:59:45W and the
last column from 179:59:45E to 180:00:15E spatially overlap, the first and
last column of DMSP are duplicates with regard to their location. If you
want to avoid this warning, you can set the region to w=179:59:45W
e=180:00:15E.
>>>
>>> Markus M
>>>
>>> >
>>> > Another question I wanted to ask is: how to know whether the
operation of intercallibration was correctly done, for tha I thought maybe
thare is the a place from where I can corroborate whether the min and max
values of each intercallibrated raster layer is correct?
>>> >
>>> >
>>> > I'm attaching the log of all the files I got from 'r.info' command in
it there appears always for the region '360 degree EW extent is exceeded by
1 cells' and  also the min and max value of each intercallibrated raster
 layer.
>>> >
>>> > So as to know if I got all the raster correctly intercallibrated
maybe checking if the min and max value for each intercallibrated
corresponds correctly is there a place where I can check that?
>>> >
>>> > Maybe according to my attached log file is possible to know if all
the intercallibration operation was correctly done and thus the layers are
ready for further study and analysis.
>>> >
>>> >
>>> > Thanks  a lot again for your help.
>>> > Kind regards,
>>> > Gabriel
>>> >
>>> > Virus-free. www.avast.com
>>> >
>>> > On Mon, Aug 19, 2019 at 4:41 PM Markus Metz <
markus.metz.giswork at gmail.com> wrote:
>>> >>
>>> >>
>>> >>
>>> >> On Mon, Aug 19, 2019 at 12:05 AM Gabriel Cotlier <gabiklm01 at gmail.com>
wrote:
>>> >> >
>>> >> > Hello,
>>> >> > My question is how does it influence the fact that it say:
>>> >> > 360 degree EW extent is exceeded by 0.999827 cells
>>> >>
>>> >> this is caused by the truncated resolution of 0.008333333300000
>>> >> with a corrected resolution of 00:00:30, the message  is
>>> >>
>>> >> > 360 degree EW extent is exceeded by 1 cells
>>> >>
>>> >> considering the EW extents of 180:00:15W to 180:00:15E, that means
that the first column from 180:00:15W to 179:59:45W and the last column
from 179:59:45E to 180:00:15E spatially overlap, the first and last column
of DMSP are duplicates with regard to their location. If you want to avoid
this warning, you can set the region to w=179:59:45W e=180:00:15E.
>>> >>
>>> >> Note that the recommended way to set a computational region to a
raster map is g.region rast=name_of_raster_map. After that, as for DMSP,
you might want to adjust the computational region to your needs, e.g. a
smaller region of interest, or restrict it to 360 degrees EW extent in case
the raster map is exceeding 360 degrees EW extent.
>>> >>
>>> >> HTH,
>>> >>
>>> >> Markus M
>>> >>
>>> >> >
>>> >> >
>>> >> >
>>> >> >
>>> >> >
>>> >> >
>>> >> > while when I loaded a first file I defined a region as
>>> >> >
>>> >> >
>>> >> > which is exactly I suppose the correct region for the DMSP
data.... then after loading the  other layers it appears:
>>> >> >
>>> >> > 360 degree EW extent is exceeded by 0.999827 cells
>>> >> > 360 degree EW extent is exceeded by 1 cells
>>> >> >
>>> >> > Thanks a lot
>>> >> > Gabriel
>>> >> >
>>> >> >
>>> >> >
>>> >> >
>>> >> > On Sun, Aug 18, 2019 at 6:54 PM Gabriel Cotlier <
gabiklm01 at gmail.com> wrote:
>>> >> >>
>>> >> >> Hello, another question, regarding
 i.nightlights.intercalibration, can I run this code as python
package/lbrary loading it from Spyder or Jupiter Notebook instead of using
GRASS interface, if so how is a convenient  way to install
i.nightlights.intercalibration in python using Spyder?
>>> >> >> Thanks a lot.
>>> >> >> Gabriel
>>> >> >>
>>> >> >> On Sat, Aug 17, 2019 at 4:54 PM Gabriel Cotlier <
gabiklm01 at gmail.com> wrote:
>>> >> >>>
>>> >> >>> Dear Nikos.
>>> >> >>> After a long time I'm trying to reproduce a routine I have for
doing intercallibratrion of DMSP 1992-2012 but for some reason It doesn't
work to me. I think is because the problem between the region of the layers
30 arc sec should resolution be from 0.008333333300000 to
0.008333333333333, i.e. exactly 30 arc-seconds? and the computational
region be the same ? I got stuck on how to set it to work... from the side
of the region setting.
>>> >> >>> However in addition my routing also has a for loop which does
not work ok as well.
>>> >> >>> I would appreciate a lot of you can give it a look and tell me
how to make it work...
>>> >> >>> Thanks  a lot in advance
>>> >> >>> Kind regards,
>>> >> >>> Gabriel
>>> >> >>>
#####-----------------------------------------------------------------------------------------
>>> >> >>> # complete routine for intercalliration of DSMP/OLS light stable
product
>>> >> >>>
>>> >> >>> import grass.script as gscript
>>> >> >>> import os
>>> >> >>> import os,glob
>>> >> >>>
>>> >> >>> # get working directory
>>> >> >>> print os.getcwd()
>>> >> >>>
>>> >> >>> # change working directory where raster files are
>>> >> >>> os.chdir('C:\\Users\\Gabriel\\Documents\\grassdata\\lights')
>>> >> >>>
>>> >> >>> # see files in directory
>>> >> >>> ls
>>> >> >>>
>>> >> >>> # import all raster files to grass --- here is a kind of
problem...???
>>> >> >>> for tif_file in glob.glob("*.tif"):
>>> >> >>>     new_rast = os.path.splitext(tif_file)[0]
>>> >> >>>     grass.run_command("r.in.gdal", flags="a", input=tif_file,
output=new_rast)
>>> >> >>>
>>> >> >>> # get info of one of the imported raster
>>> >> >>> r.info map=F121996
>>> >> >>>
>>> >> >>> # run intercalliration algorithm
>>> >> >>> i.nightlights.intercalibration
image=F101992,F101993,F101994,F121994,F121995,F121996,F121997,F121998,F121999,F141997,F141998,F141999,F142000,F142001,F142002,F142003,F152000,F152001,F152002,F152003,F152004,F152005,F152006,F152007,F162004,F162005,F162006,F162007,F162008,F162009,F182010,F182011,F182012,F182013
suffix=c model=elvidge2014 -t
>>> >> >>>
>>> >> >>> # correct general region adjust to raster file --- here the
region is exactly 30 arc for the raster as I could see....
>>> >> >>> g.region raster=F121996
>>> >> >>>
>>> >> >>> # cerate a list of rasters in the mapset
>>> >> >>> # rastlist=grass.read_command("g.list",type="rast").split()
>>> >> >>>  rasters = grass.read_command('g.list',
type='raster').splitlines()
>>> >> >>>
>>> >> >>> # change working directory
>>> >> >>> os.chdir('C:\\Users\\Gabriel\\Desktop\\out')
>>> >> >>>
>>> >> >>> # save rasters in mapset to file
>>> >> >>> for raster in rasters:
>>> >> >>>     grass.run_command('r.out.gdal', input=raster, output=raster
+ '.tiff', format='GTiff')
>>> >> >>>
>>> >> >>> On Wed, Aug 22, 2018 at 10:06 AM Gabriel Cotlier <
gabiklm01 at gmail.com> wrote:
>>> >> >>>>
>>> >> >>>> Dear Nikos,
>>> >> >>>>
>>> >> >>>> Thanks a lot for your answer and the orientation.
>>> >> >>>> The information and the link are very useful.
>>> >> >>>> Kind regards,
>>> >> >>>> Gabriel
>>> >> >>>>
>>> >> >>>>
>>> >> >>>> On Wed, Aug 22, 2018 at 5:19 AM Nikos Alexandris <
nik at nikosalexandris.net> wrote:
>>> >> >>>>>
>>> >> >>>>> * Gabriel Cotlier <gabiklm01 at gmail.com> [2018-08-21 12:00:24
-0300]:
>>> >> >>>>>
>>> >> >>>>> >Dear Nikos and GRASS users,
>>> >> >>>>> >
>>> >> >>>>> >I would like to ask if nonetheless the effect due to "stray
light" the
>>> >> >>>>> >*i.landsat8.swlst* code for split window is still applicable
to Landsat 8
>>> >> >>>>> >data and whether these error is specially visible on water
bodies? and
>>> >> >>>>> >whether band 10 is better than band 11 in terms of correction
processing
>>> >> >>>>> >for Level -1  data products?
>>> >> >>>>> >
>>> >> >>>>> >Thanks a lot.
>>> >> >>>>> >
>>> >> >>>>> >Kind regards,
>>> >> >>>>> >Gabriel
>>> >> >>>>>
>>> >> >>>>> Dear Gabriel,
>>> >> >>>>>
>>> >> >>>>> for details and references, refer to
>>> >> >>>>>
>>> >> >>>>>
https://landsat.gsfc.nasa.gov/landsat-8-thermal-data-ghost-free-after-stray-light-exorcism/
>>> >> >>>>>
>>> >> >>>>> Make sure you use the newest Level-1 Collection 1 Landsat 8
products.
>>> >> >>>>>
>>> >> >>>>> I use `i.landsat8.swlst` and plan to improve it further.
>>> >> >>>>>
>>> >> >>>>> However, whether to prefer a Split-Window based approach, or
another
>>> >> >>>>> Single-Channel one, depends on what you want to do. Think of
spatial
>>> >> >>>>> extent and coverage of various land (cover) types, temporal
extent
>>> >> >>>>> and more.
>>> >> >>>>>
>>> >> >>>>> Thermal remote sensing is hard(er) also because it's hard to
get
>>> >> >>>>> ground-truth data sets so as to validate LST estimations.
>>> >> >>>>>
>>> >> >>>>> Nikos
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-user/attachments/20190820/aacb56ec/attachment-0001.html>


More information about the grass-user mailing list