Add Physics Based Rendering (PBR) shading method to grdimage and grdview

doc/rst/source/explain_intense.rst_

@@ -14,3 +14,40 @@
     modifiers [Default is no illumination]. **Note**: If the input grid
     actually is an image then *file* or constant *intens* must be provided
     since derivatives are not available.
+
+**Or**
+
+**-I**\ **+P**\ [*azim*\ [/*elev*]][**+f**\ *fill*][**+i**\ *ior*][**+l**\ *light*][**+m**\ *metallic*][**+o**\ [*radius*]][**+r**\ *roughness*][**+s**][**+t**][**+v**\ *ve*]
+    Instead of intensities, give **+P**\ [*azim*\ [/*elev*]][*modifiers*] to shade the input grid with
+    physically based lighting (the two cannot be combined). Each node's CPT color is used as the albedo
+    of a microfacet material that is lit by the sun, at azimuth *azim* (degrees clockwise from north,
+    taken modulo 360) [315] and elevation *elev* (degrees above the horizon, 0-90) [45], and by a
+    headlight fill, and seen from straight above. Give just *azim* to keep the default *elev* (but *elev*
+    cannot be given without *azim*). Plain **+P** selects all the default settings, which are the
+    same as **+P**\ 315/45+o+t.
+
+    This is physically based rendering (PBR) as done by the shader of the Visualization Toolkit (VTK):
+    the Cook-Torrance specular term with the GGX distribution (Walter et al., 2007), the
+    height-correlated Smith visibility (Heitz, 2014) and Schlick's Fresnel approximation (Schlick,
+    1994), plus the Lambertian diffuse term that the specular leaves over. Unlike intensities, which
+    blend each color toward white or gray, the light here multiplies the color, so hues are kept in
+    both the lit and the shaded parts. The relief is lit as a 3-D view would draw it: *x* and *y* at
+    the same scale (for a geographic grid, *x* is shortened by the cosine of the middle latitude) and
+    the *z* range drawn as 0.1 of the larger horizontal side, times the vertical exaggeration (see
+    **+v**). Optionally, append modifiers (the allowed range of each value is given in parentheses,
+    the default in brackets):
+
+    - **+f** Set the *fill* (headlight) intensity (>= 0) [0.35].
+    - **+i** Set the index of refraction *ior* of the surface (>= 1) [1.5]. Values from about 1.3
+      (water) to 2.4 (diamond) cover most materials.
+    - **+l** Set the *light* (sun) intensity (>= 0) [1].
+    - **+m** Set the *metallic* value, from 0 (dielectric) to 1 (metal) (0-1) [0].
+    - **+o** Add ambient occlusion, which darkens creases and valley floors [no occlusion]. Optionally
+      append the sampling *radius* as a fraction of the diagonal of the relief as drawn (> 0) [0.1].
+      Occlusion is only computed when the radius is larger than one grid node spacing.
+    - **+r** Set the *roughness*, from 0 (mirror) to 1 (matte) (0-1) [0.3]. Values below 0.05 are
+      treated as 0.05.
+    - **+s** Cast shadows: the terrain shades itself along the sun's azimuth and elevation, which the
+      hillshade methods of :doc:`grdgradient` do not do [no shadows].
+    - **+t** Apply the Khronos PBR Neutral tone mapping, which compresses the highlights [no tone mapping].
+    - **+v** Set the vertical exaggeration *ve* of the relief (> 0) [1].

doc/rst/source/grd2kml.rst

@@ -19,7 +19,7 @@ Synopsis
 [ |-E|\ *URL* ]
 [ |-F|\ *filtercode* ]
 [ |-H|\ *scale* ]
-[ |-I|\ [*file*\|\ *intens*\|\ **+a**\ *azimuth*][**+d**][**+m**\ *ambient*][**+n**\ *args*] ]
+[ |-I|\ [*file*\|\ *intens*\|\ **+a**\ *azimuth*][**+d**][**+m**\ *ambient*][**+n**\ *args*]\|\ **+P**\ [*modifiers*] ]
 [ |-L|\ *tilesize* ]
 [ |-S|\ [*extra*] ]
 [ |-T|\ *title* ]

doc/rst/source/grdimage.rst

@@ -19,7 +19,7 @@ Synopsis
 [ |-D|\ [**r**] ]
 [ |-E|\ [**i**\|\ *dpi*] ]
 [ |-G|\ *color*\ [**+b**\|\ **f**] ]
-[ |-I|\ [*file*\|\ *intens*\|\ **+a**\ *azimuth*][**+d**][**+m**\ *ambient*][**+n**\ *args*] ]
+[ |-I|\ [*file*\|\ *intens*\|\ **+a**\ *azimuth*][**+d**][**+m**\ *ambient*][**+n**\ *args*]\|\ **+P**\ [*modifiers*] ]
 [ |-J|\ *parameters* ]
 [ |-M| ]
 [ |-N| ]
@@ -333,8 +333,26 @@ To create a sinusoidal projection of a remotely located Jessica Rabbit::
 
     gmt grdimage -JI15c -Rd http://larryfire.files.wordpress.com/2009/07/untooned_jessicarabbit.jpg -pdf jess
 
+To shade the Island of Hawaii and the surrounding sea floor with physically based lighting, with a low sun
+from the northwest that makes Mauna Kea and Mauna Loa cast their shadows, ambient occlusion, tone mapping,
+a threefold vertical exaggeration and a brighter sun and fill light than the defaults::
+
+    gmt grdimage @earth_gebco_15s -R-156.2/-154.7/18.8/20.35 -Cgeo -I+P315/35+o+s+t+v3+l1.4+f0.5 -JM15c -B -png hawaii
+
 .. include:: cpt_notes.rst_
 
+References
+----------
+
+Heitz, E., Understanding the Masking-Shadowing Function in Microfacet-Based BRDFs, Journal of
+Computer Graphics Techniques, Vol. 3, No. 2, 2014, pp. 48-107.
+
+Schlick, C., An Inexpensive BRDF Model for Physically-based Rendering, Computer Graphics Forum,
+Vol. 13, No. 3, 1994, pp. 233-246.
+
+Walter, B., Marschner, S. R., Li, H., and Torrance, K. E., Microfacet Models for Refraction through
+Rough Surfaces, Proceedings of the Eurographics Symposium on Rendering, 2007, pp. 195-206.
+
 See Also
 --------
 

doc/rst/source/grdview.rst

@@ -16,7 +16,7 @@ Synopsis
 [ |SYN_OPT-B| ]
 [ |-C|\ [*section*/]\ *master*\|\ *cpt*\|\ *color*\ :math:`_1`,\ *color*\ :math:`_2`\ [,\ *color*\ :math:`_3`\ ,...]\ [**+h**\ [*hinge*]][**+i**\ *dz*][**+u**\|\ **U**\ *unit*][**+s**\ *fname*] ]
 [ |-G|\ *drapegrid*\|\ *drapeimage* ]
-[ |-I|\ [*file*\|\ *intens*\|\ **+a**\ *azimuth*][**+d**][**+m**\ *ambient*][**+n**\ *args*] ]
+[ |-I|\ [*file*\|\ *intens*\|\ **+a**\ *azimuth*][**+d**][**+m**\ *ambient*][**+n**\ *args*]\|\ **+P**\ [*modifiers*] ]
 [ |-Jz|\ \|\ **Z**\ *parameters* ]
 [ |-N|\ [*level*]\ [**+g**\ *fill*] ]
 [ |-Q|\ **c**\|\ **g**\ [**m**]\ **i**\|\ **m**\ [**x**\|\ **y**]\|\ **s**\ [**m**]\|\ [*color*][**+m**]\ **t** ]
@@ -243,6 +243,13 @@ topo_intens.nc, and looking from the SE, run::
 
     gmt grdview topography.nc -JM6i -Gmagnetics.nc -Cmag_intens.cpt -Qs -p140/30 -Itopo_intens.nc -pdf draped3D
 
+To view the Island of Hawaii and the surrounding sea floor in perspective, shaded with physically based
+lighting: a low sun from the northwest that makes Mauna Kea and Mauna Loa cast their shadows, ambient
+occlusion, tone mapping, a threefold vertical exaggeration of the shading and a brighter sun and fill
+light than the defaults::
+
+    gmt grdview @earth_gebco_15s -R-156.2/-154.7/18.8/20.35 -Cgeo -I+P315/35+o+s+t+v3+l1.4+f0.5 -JM15c -JZ3c -p150/35 -B -png hawaii
+
 .. module_note_begins
 
 Notes
@@ -261,6 +268,18 @@ variations within polygons by using scanline conversion to image the polygons.
 
 .. module_note_ends
 
+References
+----------
+
+Heitz, E., Understanding the Masking-Shadowing Function in Microfacet-Based BRDFs, Journal of
+Computer Graphics Techniques, Vol. 3, No. 2, 2014, pp. 48-107.
+
+Schlick, C., An Inexpensive BRDF Model for Physically-based Rendering, Computer Graphics Forum,
+Vol. 13, No. 3, 1994, pp. 233-246.
+
+Walter, B., Marschner, S. R., Li, H., and Torrance, K. E., Microfacet Models for Refraction through
+Rough Surfaces, Proceedings of the Eurographics Symposium on Rendering, 2007, pp. 195-206.
+
 See Also
 --------
 

src/gmt_prototypes.h

@@ -747,6 +747,11 @@ EXTERN_MSC void gmt_RI_prepare (struct GMT_CTRL *GMT, struct GMT_GRID_HEADER *h)
 EXTERN_MSC struct GMT_DATASEGMENT * gmt_prepare_contour (struct GMT_CTRL *GMT, double *x, double *y, uint64_t n, double z);
 EXTERN_MSC void gmt_get_plot_array (struct GMT_CTRL *GMT);
 EXTERN_MSC void gmt_illuminate (struct GMT_CTRL *GMT, double intensity, double *rgb);
+EXTERN_MSC void gmt_pbr_defaults(struct GMT_PBR *P);
+EXTERN_MSC unsigned int gmt_pbr_parse(struct GMT_CTRL *GMT, char option, char *arg, struct GMT_PBR *P);
+EXTERN_MSC unsigned int gmt_pbr_check(struct GMT_CTRL *GMT, char option, struct GMT_PBR *P);
+EXTERN_MSC void gmt_pbr_syntax(struct GMTAPI_CTRL *API, char option);
+EXTERN_MSC struct GMT_IMAGE *gmt_pbr_image(struct GMT_CTRL *GMT, struct GMT_PBR *P, struct GMT_GRID *G, struct GMT_PALETTE *CPT);
 EXTERN_MSC int gmt_colorname2index (struct GMT_CTRL *GMT, char *name);
 EXTERN_MSC void gmt_list_custom_symbols (struct GMT_CTRL *GMT);
 EXTERN_MSC void gmt_smart_justify (struct GMT_CTRL *GMT, int just, double angle, double dx, double dy, double *x_shift, double *y_shift, unsigned int mode);