View source for Earth at night

== Import NASA's Night on Earth and Reproject in Polar Sterographic ==

''This micro-tutorial requires GRASS 6.5''

=== Get data ===

Download
http://veimages.gsfc.nasa.gov//1438/land_ocean_ice_lights_8192.tif
from
http://visibleearth.nasa.gov/view_detail.php?id=1438

=== Import ===
 # create simple XY location
 r.in.gdal in=land_ocean_ice_lights_8192.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
 
 #FIXME
 #hack PROJ WIND and cellhd to be 3 not 0
 #run g.setproj to set datum == wgs84 and create PROJ_INFO etc.
 # now this location is lat/lon
 
 g.region n=90 s=45 w=-180 e=180
 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 ===

 #### NORTHERN HEMISPHERE #################
 #create new loc'n 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  <>
 
 # pull in grids
 v.proj in=grid_15 location=xy_to_ll
 v.proj in=grid_30 location=xy_to_ll
 
 d.grid -w 30 color=red
 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=xy_to_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
 
 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
 
 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 ===

 #### SOUTHERN HEMISPHERE #################
 #create new loc'n based on inverse of epsg #3413
 #proj       : stere
 #datum      : wgs84
 #lat_0      : 90
 #lat_ts     : 70
 #lon_0      : +135
 #k          : 1
 
 v.proj in=grid_15 location=xy_to_ll
 v.proj in=grid_30 location=xy_to_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=xy_to_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
 
 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
 
 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


=== Result ===

* {{website|screenshots/images/night_on_earth_N.png|Northern Hemisphere}}
* {{website|screenshots/images/night_on_earth_S.png|Southern Hemisphere}}


[[Category: Documentation]]
[[Category: Tutorial]]