[GRASS-dev] [GRASS GIS] #2014: r.sun using EPSG:3031 projection gives strange results

[GRASS GIS]

#2014: r.sun using EPSG:3031 projection gives strange results
---------------------------+------------------------------------------------
 Reporter:  pierreroudier  |       Owner:  grass-dev@…              
     Type:  defect         |      Status:  new                      
 Priority:  normal         |   Milestone:  7.0.0                    
Component:  Raster         |     Version:  svn-trunk                
 Keywords:  r.sun          |    Platform:  Linux                    
      Cpu:  x86-64         |  
---------------------------+------------------------------------------------

Comment(by pierreroudier):

 Hamish,

 The original command I used was:

 {{{

 ## annual beam radiance
 DEM=rsr_dem

 # setup region, and re-compute slope in deg.
  g.region rast=$DEM res=200 -ap
 #g.region res=200 n=-567200 s=-634800 w=-3400 e=117600 -ap

 r.slope.aspect --o elev=$DEM slope=slope_deg aspect=aspect

 # setup: no accomodation for Linke, using defaults
 BEGIN=1
 END=365
 STEP=1
 NUM_CORES=4

 for DAY in `seq $BEGIN $STEP $END` ; do
 DAY_STR=`echo $DAY | awk '{printf("%.03d", $1)}'`
 echo "Processing day $DAY_STR at `date` ..."
 echo "Processing day $DAY_STR at `date` ...\n" >> r.sun_progress.txt

 CMD="r.sun --o --quiet \
          elev_in=$DEM slope_in=slope_deg asp_in=aspect \
          step=1.0 \
          day=$DAY \
          beam_rad=beam_rad.$DAY_STR \
          insol_time=insol_time.$DAY_STR "

 # poor man's multi-threading for a multi-core CPU
 MODULUS=`echo "$DAY $STEP $NUM_CORES" | awk '{print $1 % ($2 * $3)}'`
 if [ "$MODULUS" = "$STEP" ] ; then
    # stall to let the background jobs finish
    $CMD
    sleep 2
 else
    $CMD &
 fi
 done


 # combine into annual sum
 r.series --o in=`g.mlist type=rast pattern="beam_rad.*" sep=","`
 out=beam_rad_annual method=sum
 r.series --o in=`g.mlist type=rast pattern="insol_time.*" sep=","`
 out=insol_time_annual method=sum

 # convert to MJ/sq. m
 # r.sun mode 2 * (3600 / 100000) ==> MJ/sq. m
 r.mapcalc --o "beam_rad_annual_mj = beam_rad_annual * (3600.0 / 100000.0)"
 r.colors.stddev beam_rad_annual_mj


 # clean-up
 # g.mremove -f rast=beam_rad.*
 # g.mremove -f rast=insol_time.*
 }}}

 Happy to assist with the testing on grass7. Which Antarctic projection
 would you suggests? I got no projection in mind that would suit our needs.

 Pierre

 BTW the screenshot was generated in QGIS, not GRASS, so that's one less
 bug to fix :)

-- 
Ticket URL: <https://trac.osgeo.org/grass/ticket/2014#comment:7>
GRASS GIS <http://grass.osgeo.org>