converting the projection

Previous Topic Next Topic
 
classic Classic list List threaded Threaded
8 messages Options
Reply | Threaded
Open this post in threaded view
|

converting the projection

Kyongho Son

Hi I have two raster data sets but they have two datum and projection.

I would like to convert one to the other since I would like to clip the one raster using the other raster.

they should have same projections.


proj4string(bb)
[1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
proj4string(SOSTavg)

 proj4string(SOSTavg)
[1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"


Thanks,

Kyongho








        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting the projection

Michael Sumner-2
Try projecting just the extent
library(raster)
ex <- projectExtent(a,  projection(b))
clip(b,  ex)

where a,  b are your two raster.

Cheers,  Mike

On Wed, 20 Sep 2017, 06:18 Kyongho Son <[hidden email]> wrote:

>
> Hi I have two raster data sets but they have two datum and projection.
>
> I would like to convert one to the other since I would like to clip the
> one raster using the other raster.
>
> they should have same projections.
>
>
> proj4string(bb)
> [1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
> proj4string(SOSTavg)
>
>  proj4string(SOSTavg)
> [1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +ellps=WGS84
> +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"
>
>
> Thanks,
>
> Kyongho
>
>
>
>
>
>
>
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
--
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting the projection

Michael Sumner-2
Keep on list please.

Share code to demonstrate so we can see the problem.  Use projectRaster(a,
b)  for alignment,  try crop(b,  ex,  snap = "out")  for inclusive crop.

And apologies I wrote clip instead of crop earlier.

Cheers,  Mike

On Wed, 20 Sep 2017, 07:42 Kyongho Son <[hidden email]> wrote:

> I did the exact same way, but the output of the reprojection of 'a' raster
> did not properly align with 'b' raster.
>
>
>
>
>
> *Kyongho Son*
> Postdoctoral Fellow
> Research Foundation of The City University of New York
>
>
>
> ------------------------------
> *From:* Michael Sumner <[hidden email]>
> *Sent:* Tuesday, September 19, 2017 5:03 PM
> *To:* Kyongho Son; r-sig-geo
> *Subject:* Re: [R-sig-Geo] converting the projection
>
> Try projecting just the extent
> library(raster)
> ex <- projectExtent(a,  projection(b))
> clip(b,  ex)
>
> where a,  b are your two raster.
>
> Cheers,  Mike
>
> On Wed, 20 Sep 2017, 06:18 Kyongho Son <[hidden email]> wrote:
>
>>
>> Hi I have two raster data sets but they have two datum and projection.
>>
>> I would like to convert one to the other since I would like to clip the
>> one raster using the other raster.
>>
>> they should have same projections.
>>
>>
>> proj4string(bb)
>> [1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
>> proj4string(SOSTavg)
>>
>>  proj4string(SOSTavg)
>> [1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +ellps=WGS84
>> +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"
>>
>>
>> Thanks,
>>
>> Kyongho
>>
>>
>>
>>
>>
>>
>>
>>
>>         [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
> --
>
> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
> Kingston Tasmania 7050 Australia
> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>
> --
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting the projection

Kyongho Son
Sorry and thanks, Mike
I can not share all codes and data since the original data are huge.
I copied the scripts that I used to reproject one of rasters and then make sure that they have same projection

 Then clip a large raster with a small raster. However, I got the error messages, the two raster maps are not overlap.

They should overlap if the projection was correct.



# east US map tree phenology map
tmp<-raster("../../obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114")

proj4string(tmp)
[1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

# small watershed boundary map in NY state
bb<-raster("../../obs/gis/bb_soil_def_map.tif")
proj4string(bb)
[1] "+proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

ex <- projectExtent(bb,  projection(tmp))

> ex
class       : RasterLayer
dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
resolution  : 30.12012, 14.94084  (x, y)
extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0


c<-crop(tmp,  ex)

Error in .local(x, y, ...) : extents do not overlap




*Kyongho Son*
Postdoctoral Fellow
Research Foundation of The City University of New York


________________________________
From: Michael Sumner <[hidden email]>
Sent: Tuesday, September 19, 2017 8:32 PM
To: Kyongho Son; RsigGeo
Subject: Re: [R-sig-Geo] converting the projection

Keep on list please.

Share code to demonstrate so we can see the problem.  Use projectRaster(a, b)  for alignment,  try crop(b,  ex,  snap = "out")  for inclusive crop.

And apologies I wrote clip instead of crop earlier.

Cheers,  Mike

On Wed, 20 Sep 2017, 07:42 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:

I did the exact same way, but the output of the reprojection of 'a' raster did not properly align with 'b' raster.





*Kyongho Son*
Postdoctoral Fellow
Research Foundation of The City University of New York


________________________________
From: Michael Sumner <[hidden email]<mailto:[hidden email]>>
Sent: Tuesday, September 19, 2017 5:03 PM
To: Kyongho Son; r-sig-geo
Subject: Re: [R-sig-Geo] converting the projection

Try projecting just the extent
library(raster)
ex <- projectExtent(a,  projection(b))
clip(b,  ex)

where a,  b are your two raster.

Cheers,  Mike

On Wed, 20 Sep 2017, 06:18 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:

Hi I have two raster data sets but they have two datum and projection.

I would like to convert one to the other since I would like to clip the one raster using the other raster.

they should have same projections.


proj4string(bb)
[1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
proj4string(SOSTavg)

 proj4string(SOSTavg)
[1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"


Thanks,

Kyongho








        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]<mailto:[hidden email]>
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
Kingston Tasmania 7050 Australia<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>

--
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting the projection

Michael Sumner-2
Great thanks, that is good detail -but we still need the extents of both
data to be able to see where the problem is.

You gave us the proj4string for both, but we also need the extent - you can
print each object for a full summary eg.

print(bb)

print(tmp)

and with those anyone can reproduce enough (without any real data) to
recreate your experience. (I am so used to working with this stuff I forget
to remind others that print(theraster) is usually more than enough
information to work through a reproducible example)

But, a guess is - unless the projection metadata is wrong, raster is
probably telling you the truth - but it should be easy to see why there's a
problem with the print output (essentially the proj4string and extent
provide enough information for an experienced problem solver).

Cheers, Mike.



On Wed, 20 Sep 2017 at 11:05 Kyongho Son <[hidden email]> wrote:

> Sorry and thanks, Mike
> I can not share all codes and data since the original data are huge.
> I copied the scripts that I used to reproject one of rasters and then make
> sure that they have same projection
>
>  Then clip a large raster with a small raster. However, I got the error
> messages, the two raster maps are not overlap.
>
> They should overlap if the projection was correct.
>
>
>
> # east US map tree phenology map
> tmp<-raster("../../obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114")
>
> proj4string(tmp)
> [1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84
> +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
>
> # small watershed boundary map in NY state
> bb<-raster("../../obs/gis/bb_soil_def_map.tif")
> proj4string(bb)
> [1] "+proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84
> +towgs84=0,0,0"
>
> ex <- projectExtent(bb,  projection(tmp))
>
> > ex
> class       : RasterLayer
> dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
> resolution  : 30.12012, 14.94084  (x, y)
> extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
> coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0
> +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
>
>
> c<-crop(tmp,  ex)
>
> Error in .local(x, y, ...) : extents do not overlap
>
>
>
> *Kyongho Son*
> Postdoctoral Fellow
> Research Foundation of The City University of New York
>
>
>
> ------------------------------
> *From:* Michael Sumner <[hidden email]>
> *Sent:* Tuesday, September 19, 2017 8:32 PM
> *To:* Kyongho Son; RsigGeo
>
> *Subject:* Re: [R-sig-Geo] converting the projection
> Keep on list please.
>
> Share code to demonstrate so we can see the problem.  Use projectRaster(a,
> b)  for alignment,  try crop(b,  ex,  snap = "out")  for inclusive crop.
>
> And apologies I wrote clip instead of crop earlier.
>
> Cheers,  Mike
>
> On Wed, 20 Sep 2017, 07:42 Kyongho Son <[hidden email]> wrote:
>
>> I did the exact same way, but the output of the reprojection of 'a'
>> raster did not properly align with 'b' raster.
>>
>>
>>
>>
>>
>> *Kyongho Son*
>> Postdoctoral Fellow
>> Research Foundation of The City University of New York
>>
>>
>>
>> ------------------------------
>> *From:* Michael Sumner <[hidden email]>
>> *Sent:* Tuesday, September 19, 2017 5:03 PM
>> *To:* Kyongho Son; r-sig-geo
>> *Subject:* Re: [R-sig-Geo] converting the projection
>>
>> Try projecting just the extent
>> library(raster)
>> ex <- projectExtent(a,  projection(b))
>> clip(b,  ex)
>>
>> where a,  b are your two raster.
>>
>> Cheers,  Mike
>>
>> On Wed, 20 Sep 2017, 06:18 Kyongho Son <[hidden email]> wrote:
>>
>>>
>>> Hi I have two raster data sets but they have two datum and projection.
>>>
>>> I would like to convert one to the other since I would like to clip the
>>> one raster using the other raster.
>>>
>>> they should have same projections.
>>>
>>>
>>> proj4string(bb)
>>> [1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
>>> proj4string(SOSTavg)
>>>
>>>  proj4string(SOSTavg)
>>> [1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +ellps=WGS84
>>> +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"
>>>
>>>
>>> Thanks,
>>>
>>> Kyongho
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>         [[alternative HTML version deleted]]
>>>
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> [hidden email]
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>> --
>>
>> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g>
>> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>> Dr. Michael Sumner
>> Software and Database Engineer
>> Australian Antarctic Division
>> 203 Channel Highway
>> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>> Kingston Tasmania 7050 Australia
>> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>>
>> --
>
> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g>
> Kingston Tasmania 7050 Australia
> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g>
>
> --
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting the projection

Kyongho Son
I add the extent

# reading emodis data

# east US map tree phenology map
tmp<-raster("../../obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114")

proj4string(tmp)
[1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

print(tmp)

> print(tmp)
class       : RasterLayer
dimensions  : 11556, 9174, 106014744  (nrow, ncol, ncell)
resolution  : 250, 250  (x, y)
extent      : 243000, 2536500, -2136500, 752500  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
data source : /work.taichi/Projects/Biscuit/rhessys/obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114
names       : sostavg0114
values      : -115, 1000  (min, max)


# small watershed boundary map in NY state
bb<-raster("../../obs/gis/bb_soil_def_map.tif")
proj4string(bb)
[1] "+proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

print(bb)

> print(bb)
class       : RasterLayer
dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
resolution  : 10, 10  (x, y)
extent      : 540217, 545367, 4648219, 4654199  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
data source : /work.taichi/Projects/Biscuit/rhessys/obs/gis/bb_soil_def_map.tif
names       : bb_soil_def_map
values      : 0, 5712  (min, max)

ex <- projectExtent(bb,  projection(tmp))

print(ex)

> print(ex)
class       : RasterLayer
dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
resolution  : 30.12012, 14.94084  (x, y)
extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0

> ex
class       : RasterLayer
dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
resolution  : 30.12012, 14.94084  (x, y)
extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0


c<-crop(tmp,  ex)

Error in .local(x, y, ...) : extents do not overlap



________________________________
From: Michael Sumner <[hidden email]>
Sent: Tuesday, September 19, 2017 9:15 PM
To: Kyongho Son; RsigGeo
Subject: Re: [R-sig-Geo] converting the projection

Great thanks, that is good detail -but we still need the extents of both data to be able to see where the problem is.

You gave us the proj4string for both, but we also need the extent - you can print each object for a full summary eg.

print(bb)

print(tmp)

and with those anyone can reproduce enough (without any real data) to recreate your experience. (I am so used to working with this stuff I forget to remind others that print(theraster) is usually more than enough information to work through a reproducible example)

But, a guess is - unless the projection metadata is wrong, raster is probably telling you the truth - but it should be easy to see why there's a problem with the print output (essentially the proj4string and extent provide enough information for an experienced problem solver).

Cheers, Mike.



On Wed, 20 Sep 2017 at 11:05 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:

Sorry and thanks, Mike
I can not share all codes and data since the original data are huge.
I copied the scripts that I used to reproject one of rasters and then make sure that they have same projection

 Then clip a large raster with a small raster. However, I got the error messages, the two raster maps are not overlap.

They should overlap if the projection was correct.



# east US map tree phenology map
tmp<-raster("../../obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114")

proj4string(tmp)
[1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

# small watershed boundary map in NY state
bb<-raster("../../obs/gis/bb_soil_def_map.tif")
proj4string(bb)
[1] "+proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

ex <- projectExtent(bb,  projection(tmp))

> ex
class       : RasterLayer
dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
resolution  : 30.12012, 14.94084  (x, y)
extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0


c<-crop(tmp,  ex)

Error in .local(x, y, ...) : extents do not overlap




*Kyongho Son*
Postdoctoral Fellow
Research Foundation of The City University of New York


________________________________
From: Michael Sumner <[hidden email]<mailto:[hidden email]>>
Sent: Tuesday, September 19, 2017 8:32 PM
To: Kyongho Son; RsigGeo

Subject: Re: [R-sig-Geo] converting the projection
Keep on list please.

Share code to demonstrate so we can see the problem.  Use projectRaster(a, b)  for alignment,  try crop(b,  ex,  snap = "out")  for inclusive crop.

And apologies I wrote clip instead of crop earlier.

Cheers,  Mike

On Wed, 20 Sep 2017, 07:42 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:

I did the exact same way, but the output of the reprojection of 'a' raster did not properly align with 'b' raster.





*Kyongho Son*
Postdoctoral Fellow
Research Foundation of The City University of New York


________________________________
From: Michael Sumner <[hidden email]<mailto:[hidden email]>>
Sent: Tuesday, September 19, 2017 5:03 PM
To: Kyongho Son; r-sig-geo
Subject: Re: [R-sig-Geo] converting the projection

Try projecting just the extent
library(raster)
ex <- projectExtent(a,  projection(b))
clip(b,  ex)

where a,  b are your two raster.

Cheers,  Mike

On Wed, 20 Sep 2017, 06:18 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:

Hi I have two raster data sets but they have two datum and projection.

I would like to convert one to the other since I would like to clip the one raster using the other raster.

they should have same projections.


proj4string(bb)
[1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
proj4string(SOSTavg)

 proj4string(SOSTavg)
[1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"


Thanks,

Kyongho








        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]<mailto:[hidden email]>
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g><https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
Kingston Tasmania 7050 Australia<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>

--
<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g>
Kingston Tasmania 7050 Australia<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g>

--
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

converting the projection from Lambert to WGS84 or NAD83

Kyongho Son

I am trying to reproject the eMODIS raster

  *   Projected Coordinate System:Sphere_ARC_INFO_Lambert_Azimuthal_Equal_Area
  *   Projection: Lambert_Azimuthal_Equal_Area
  *   Datum: D_Sphere_ARC_INFO

to the WGS84 or NAD83 18N

So that I can clip eMODIS raster files (more than 100) using a watershed boundary.
I attached the script here, but I still can not properly reproject the eMOIDS raster files.


# reading emodis data

# east US map tree phenology map
tmp<-raster("../../obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114")

proj4string(tmp)
[1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

print(tmp)

> print(tmp)
class       : RasterLayer
dimensions  : 11556, 9174, 106014744  (nrow, ncol, ncell)
resolution  : 250, 250  (x, y)
extent      : 243000, 2536500, -2136500, 752500  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
data source : /work.taichi/Projects/Biscuit/rhessys/obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114
names       : sostavg0114
values      : -115, 1000  (min, max)


# small watershed boundary map in NY state
bb<-raster("../../obs/gis/bb_soil_def_map.tif")
proj4string(bb)
[1] "+proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

print(bb)

> print(bb)
class       : RasterLayer
dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
resolution  : 10, 10  (x, y)
extent      : 540217, 545367, 4648219, 4654199  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
data source : /work.taichi/Projects/Biscuit/rhessys/obs/gis/bb_soil_def_map.tif
names       : bb_soil_def_map
values      : 0, 5712  (min, max)

ex <- projectExtent(bb,  projection(tmp))

print(ex)

> print(ex)
class       : RasterLayer
dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
resolution  : 30.12012, 14.94084  (x, y)
extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0

> ex
class       : RasterLayer
dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
resolution  : 30.12012, 14.94084  (x, y)
extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0


c<-crop(tmp,  ex)

Error in .local(x, y, ...) : extents do not overlap



________________________________
From: Michael Sumner <[hidden email]>
Sent: Tuesday, September 19, 2017 9:15 PM
To: Kyongho Son; RsigGeo
Subject: Re: [R-sig-Geo] converting the projection

Great thanks, that is good detail -but we still need the extents of both data to be able to see where the problem is.

You gave us the proj4string for both, but we also need the extent - you can print each object for a full summary eg.

print(bb)

print(tmp)

and with those anyone can reproduce enough (without any real data) to recreate your experience. (I am so used to working with this stuff I forget to remind others that print(theraster) is usually more than enough information to work through a reproducible example)

But, a guess is - unless the projection metadata is wrong, raster is probably telling you the truth - but it should be easy to see why there's a problem with the print output (essentially the proj4string and extent provide enough information for an experienced problem solver).

Cheers, Mike.



On Wed, 20 Sep 2017 at 11:05 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:

Sorry and thanks, Mike
I can not share all codes and data since the original data are huge.
I copied the scripts that I used to reproject one of rasters and then make sure that they have same projection

 Then clip a large raster with a small raster. However, I got the error messages, the two raster maps are not overlap.

They should overlap if the projection was correct.



# east US map tree phenology map
tmp<-raster("../../obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114")

proj4string(tmp)
[1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

# small watershed boundary map in NY state
bb<-raster("../../obs/gis/bb_soil_def_map.tif")
proj4string(bb)
[1] "+proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

ex <- projectExtent(bb,  projection(tmp))

> ex
class       : RasterLayer
dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
resolution  : 30.12012, 14.94084  (x, y)
extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0


c<-crop(tmp,  ex)

Error in .local(x, y, ...) : extents do not overlap




*Kyongho Son*
Postdoctoral Fellow
Research Foundation of The City University of New York


________________________________
From: Michael Sumner <[hidden email]<mailto:[hidden email]>>
Sent: Tuesday, September 19, 2017 8:32 PM
To: Kyongho Son; RsigGeo

Subject: Re: [R-sig-Geo] converting the projection
Keep on list please.

Share code to demonstrate so we can see the problem.  Use projectRaster(a, b)  for alignment,  try crop(b,  ex,  snap = "out")  for inclusive crop.

And apologies I wrote clip instead of crop earlier.

Cheers,  Mike

On Wed, 20 Sep 2017, 07:42 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:

I did the exact same way, but the output of the reprojection of 'a' raster did not properly align with 'b' raster.





*Kyongho Son*
Postdoctoral Fellow
Research Foundation of The City University of New York


________________________________
From: Michael Sumner <[hidden email]<mailto:[hidden email]>>
Sent: Tuesday, September 19, 2017 5:03 PM
To: Kyongho Son; r-sig-geo
Subject: Re: [R-sig-Geo] converting the projection

Try projecting just the extent
library(raster)
ex <- projectExtent(a,  projection(b))
clip(b,  ex)

where a,  b are your two raster.

Cheers,  Mike

On Wed, 20 Sep 2017, 06:18 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:

Hi I have two raster data sets but they have two datum and projection.

I would like to convert one to the other since I would like to clip the one raster using the other raster.

they should have same projections.


proj4string(bb)
[1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
proj4string(SOSTavg)

 proj4string(SOSTavg)
[1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"


Thanks,

Kyongho








        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]<mailto:[hidden email]>
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g><https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
Kingston Tasmania 7050 Australia<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>

--
<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g>
Kingston Tasmania 7050 Australia<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g>

--
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting the projection from Lambert to WGS84 or NAD83

Francisco Rodriguez Sanchez
Hi Kyongho,

There seems to be something strange with your eMODIS raster, which I
understood represents something in the eastern US. There might be
something wrong with the coordinates or the original projection
information. When I reprojected it to global lonlat to see which area it
covers, it actually fell in Antartica! (see
https://i.imgur.com/wjTBt8E.png). So, I'm not sure if I got it
correctly, but I think that's not what you expect. That would explain
the error that the extents do not overlap.

Reproducible example follows:

library(raster)

## Create raster with same extent and projection
emodis <- raster(xmn = 243000, xmx = 2536500, ymn = -2136500, ymx = 752500,
                  resolution = 10000,
                  crs = "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45
+y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0",
                  vals = 0)

emodis
#> class       : RasterLayer
#> dimensions  : 289, 229, 66181  (nrow, ncol, ncell)
#> resolution  : 10000, 10000  (x, y)
#> extent      : 243000, 2533000, -2137500, 752500  (xmin, xmax, ymin, ymax)
#> coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0
+datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
#> data source : in memory
#> names       : layer
#> values      : 0, 0  (min, max)

## Project to lonlat to see where it falls
emodis.lonlat <- projectRaster(emodis, crs = "+proj=longlat +datum=WGS84")

library(maptools)
data(wrld_simpl)
plot(wrld_simpl)
plot(emodis.lonlat, add = TRUE)

# See map at ![](https://i.imgur.com/wjTBt8E.png)

Hope this helps!

Cheers,

Paco



El 20/09/2017 a las 21:50, Kyongho Son escribió:

> I am trying to reproject the eMODIS raster
>
>    *   Projected Coordinate System:Sphere_ARC_INFO_Lambert_Azimuthal_Equal_Area
>    *   Projection: Lambert_Azimuthal_Equal_Area
>    *   Datum: D_Sphere_ARC_INFO
>
> to the WGS84 or NAD83 18N
>
> So that I can clip eMODIS raster files (more than 100) using a watershed boundary.
> I attached the script here, but I still can not properly reproject the eMOIDS raster files.
>
>
> # reading emodis data
>
> # east US map tree phenology map
> tmp<-raster("../../obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114")
>
> proj4string(tmp)
> [1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
>
> print(tmp)
>
>> print(tmp)
> class       : RasterLayer
> dimensions  : 11556, 9174, 106014744  (nrow, ncol, ncell)
> resolution  : 250, 250  (x, y)
> extent      : 243000, 2536500, -2136500, 752500  (xmin, xmax, ymin, ymax)
> coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
> data source : /work.taichi/Projects/Biscuit/rhessys/obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114
> names       : sostavg0114
> values      : -115, 1000  (min, max)
>
>
> # small watershed boundary map in NY state
> bb<-raster("../../obs/gis/bb_soil_def_map.tif")
> proj4string(bb)
> [1] "+proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
>
> print(bb)
>
>> print(bb)
> class       : RasterLayer
> dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
> resolution  : 10, 10  (x, y)
> extent      : 540217, 545367, 4648219, 4654199  (xmin, xmax, ymin, ymax)
> coord. ref. : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
> data source : /work.taichi/Projects/Biscuit/rhessys/obs/gis/bb_soil_def_map.tif
> names       : bb_soil_def_map
> values      : 0, 5712  (min, max)
>
> ex <- projectExtent(bb,  projection(tmp))
>
> print(ex)
>
>> print(ex)
> class       : RasterLayer
> dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
> resolution  : 30.12012, 14.94084  (x, y)
> extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
> coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
>
>> ex
> class       : RasterLayer
> dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
> resolution  : 30.12012, 14.94084  (x, y)
> extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
> coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
>
>
> c<-crop(tmp,  ex)
>
> Error in .local(x, y, ...) : extents do not overlap
>
>
>
> ________________________________
> From: Michael Sumner <[hidden email]>
> Sent: Tuesday, September 19, 2017 9:15 PM
> To: Kyongho Son; RsigGeo
> Subject: Re: [R-sig-Geo] converting the projection
>
> Great thanks, that is good detail -but we still need the extents of both data to be able to see where the problem is.
>
> You gave us the proj4string for both, but we also need the extent - you can print each object for a full summary eg.
>
> print(bb)
>
> print(tmp)
>
> and with those anyone can reproduce enough (without any real data) to recreate your experience. (I am so used to working with this stuff I forget to remind others that print(theraster) is usually more than enough information to work through a reproducible example)
>
> But, a guess is - unless the projection metadata is wrong, raster is probably telling you the truth - but it should be easy to see why there's a problem with the print output (essentially the proj4string and extent provide enough information for an experienced problem solver).
>
> Cheers, Mike.
>
>
>
> On Wed, 20 Sep 2017 at 11:05 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:
>
> Sorry and thanks, Mike
> I can not share all codes and data since the original data are huge.
> I copied the scripts that I used to reproject one of rasters and then make sure that they have same projection
>
>   Then clip a large raster with a small raster. However, I got the error messages, the two raster maps are not overlap.
>
> They should overlap if the projection was correct.
>
>
>
> # east US map tree phenology map
> tmp<-raster("../../obs/eMODIS_East_SOSTAvg0114_v1/sostavg0114")
>
> proj4string(tmp)
> [1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
>
> # small watershed boundary map in NY state
> bb<-raster("../../obs/gis/bb_soil_def_map.tif")
> proj4string(bb)
> [1] "+proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
>
> ex <- projectExtent(bb,  projection(tmp))
>
>> ex
> class       : RasterLayer
> dimensions  : 598, 515, 307970  (nrow, ncol, ncell)
> resolution  : 30.12012, 14.94084  (x, y)
> extent      : 6691756, 6707267, 9941473, 9950408  (xmin, xmax, ymin, ymax)
> coord. ref. : +proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
>
>
> c<-crop(tmp,  ex)
>
> Error in .local(x, y, ...) : extents do not overlap
>
>
>
>
> *Kyongho Son*
> Postdoctoral Fellow
> Research Foundation of The City University of New York
>
>
> ________________________________
> From: Michael Sumner <[hidden email]<mailto:[hidden email]>>
> Sent: Tuesday, September 19, 2017 8:32 PM
> To: Kyongho Son; RsigGeo
>
> Subject: Re: [R-sig-Geo] converting the projection
> Keep on list please.
>
> Share code to demonstrate so we can see the problem.  Use projectRaster(a, b)  for alignment,  try crop(b,  ex,  snap = "out")  for inclusive crop.
>
> And apologies I wrote clip instead of crop earlier.
>
> Cheers,  Mike
>
> On Wed, 20 Sep 2017, 07:42 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:
>
> I did the exact same way, but the output of the reprojection of 'a' raster did not properly align with 'b' raster.
>
>
>
>
>
> *Kyongho Son*
> Postdoctoral Fellow
> Research Foundation of The City University of New York
>
>
> ________________________________
> From: Michael Sumner <[hidden email]<mailto:[hidden email]>>
> Sent: Tuesday, September 19, 2017 5:03 PM
> To: Kyongho Son; r-sig-geo
> Subject: Re: [R-sig-Geo] converting the projection
>
> Try projecting just the extent
> library(raster)
> ex <- projectExtent(a,  projection(b))
> clip(b,  ex)
>
> where a,  b are your two raster.
>
> Cheers,  Mike
>
> On Wed, 20 Sep 2017, 06:18 Kyongho Son <[hidden email]<mailto:[hidden email]>> wrote:
>
> Hi I have two raster data sets but they have two datum and projection.
>
> I would like to convert one to the other since I would like to clip the one raster using the other raster.
>
> they should have same projections.
>
>
> proj4string(bb)
> [1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
> proj4string(SOSTavg)
>
>   proj4string(SOSTavg)
> [1] "+proj=laea +lat_0=-100 +lon_0=6370997 +x_0=45 +y_0=0 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"
>
>
> Thanks,
>
> Kyongho
>
>
>
>
>
>
>
>
>          [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]<mailto:[hidden email]>
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> --
> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g><https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
> Kingston Tasmania 7050 Australia<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>
> --
> <https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g>
> Kingston Tasmania 7050 Australia<https://maps.google.com/?q=203+Channel+Highway+%0D+Kingston+Tasmania+7050+Australia+%3Chttps://maps.google.com/?q%3D203%2BChannel%2BHighway%2B%250D%2BKingston%2BTasmania%2B7050%2BAustralia%26entry%3Dgmail%26source%3Dg%3E&entry=gmail&source=g>
>
> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> Kingston Tasmania 7050 Australia
>
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

--
Dr Francisco Rodriguez-Sanchez
Integrative Ecology Group
Estacion Biologica de Doñana (CSIC)
Avda. Americo Vespucio 26
E-41092 Sevilla (Spain)
http://bit.ly/frod_san

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo