Fractional Brownian field or surface generator

Simulates Fractional Brownian field on unit disk, with Hurst parameter 'H';
Note that the covariance function is isotropic, see reference below.
INPUTS:
- 'H' is the Hurst parameter of the Gaussian process
- 'n' is the number of grid points, where 'n' is a power of 2;
if the 'n' supplied is not a power of two,
then we set n=2^ceil(log2(n)); default is n=2^8;
OUTPUT:
- two statistically independent fields 'field1' and 'field2'
over unit disk; if not output requested, then function
outputs a figure of one of the fields
- vectors 'tx' and 'ty' so that the field is plotted via
surf(tx,ty,field1,'EdgeColor','none')
Example:
[field1,field2,tx,ty]=Brownian_field(.9,2^10);surf(tx,ty,field2,'EdgeColor','none'),colormap bone
Reference:
Kroese, D. P., & Botev, Z. I. (2015). Spatial Process Simulation.
In Stochastic Geometry, Spatial Statistics and Random Fields(pp. 369-404)
Springer International Publishing, DOI: 10.1007/978-3-319-10064-7_12