Earth at night
Jump to navigation
Jump to search
Import NASA's Night on Earth and Reproject in Polar Stereographic
This micro-tutorial requires GRASS 6.5
Get data
Download http://eoimages.gsfc.nasa.gov/ve/1438/land_ocean_ice_lights_2048.tif from NASA's Visible Earth website.
Import
- Create a WGS84 Lat/Lon location (EPSG # 4326)
# import into lat/lon, overriding projection check and # setting dummy region bounds which fit into geographic # space (the file is missing georeferencing metadata but # we do know exactly what it should be). r.in.gdal -o -l in=land_ocean_ice_lights_2048.tif out=land_ocean_ice_lights for BAND in red green blue ; do r.region map=land_ocean_ice_lights.$BAND n=90 s=-90 w=-180 e=180 done g.region n=90 s=45 w=-180 e=180 # create 15° and 30° grid line vectors v.mkgrid map=grid_15 box=15,15 grid=12,24 pos=coor coor=-180,-90 breaks=30 v.mkgrid map=grid_30 box=30,30 grid=6,12 pos=coor coor=-180,-90 breaks=30
Reproject North side
- Create new location based on EPSG #3413
# WGS 84 / NSIDC Sea Ice Polar Stereographic North <3413> +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +k=1 +x_0=0 \ +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs <>
#### NORTHERN HEMISPHERE ################# # pull in grids v.proj in=grid_15 location=world_ll v.proj in=grid_30 location=world_ll g.region n=10000000 s=-10000000 w=-10000000 e=10000000 res=5000 -p MAP=land_ocean_ice_lights for BAND in red green blue ; do echo "Processing [$BAND] ..." r.proj -n in=$MAP.$BAND location=world_ll \ method=cubic_f memory=300 r.mapcalc "${MAP}_$BAND = if($MAP.$BAND > 255, 255, if($MAP.$BAND < 0, 0, int($MAP.$BAND)))" g.remove "$MAP.$BAND" g.rename "${MAP}_$BAND,$MAP.$BAND" r.colors $MAP.$BAND color=grey255 echo done # draw d.erase d.rgb r=$MAP.red g=$MAP.green b=$MAP.blue d.vect grid_15 type=boundary color=80:40:40 d.vect grid_30 type=boundary color=120:50:50 # output d.out.file night_on_earth_N -c size=1000,1000 format=png d.out.file night_on_earth_N_small -c size=250,250 format=png
Reproject South side
- Create new location based on the inverse of EPSG #3413
+proj : stere +datum : wgs84 +lat_0 : -90 +lat_ts : -70 +lon_0 : +135 +k : 1
#### SOUTHERN HEMISPHERE ################# # pull in grids v.proj in=grid_15 location=world_ll v.proj in=grid_30 location=world_ll g.region n=10000000 s=-10000000 w=-10000000 e=10000000 res=5000 -p MAP=land_ocean_ice_lights for BAND in red green blue ; do echo "Processing [$BAND] ..." r.proj -n in=$MAP.$BAND location=world_ll \ method=cubic_f memory=600 r.mapcalc "${MAP}_$BAND = if($MAP.$BAND > 255, 255, if($MAP.$BAND < 0, 0, int($MAP.$BAND)))" g.remove "$MAP.$BAND" g.rename "${MAP}_$BAND,$MAP.$BAND" r.colors $MAP.$BAND color=grey255 echo done # remove stray NULLs for BAND in green blue ; do echo "Processing [$BAND] ..." g.rename "$MAP.$BAND,$MAP.$BAND.tmp1" r.neighbors in="$MAP.$BAND.tmp1" out="$MAP.$BAND.tmp2" r.patch in="$MAP.$BAND.tmp1,$MAP.$BAND.tmp2" out="$MAP.$BAND" g.remove "$MAP.$BAND.tmp1,$MAP.$BAND.tmp2" echo done # draw d.erase d.rgb r=$MAP.red g=$MAP.green b=$MAP.blue d.vect grid_15 type=boundary color=80:40:40 d.vect grid_30 type=boundary color=120:50:50 # output d.out.file night_on_earth_S -c size=1000,1000 format=png d.out.file night_on_earth_S_small -c size=250,250 format=png