Re: Simulating random fields in Mathematica
- To: mathgroup at smc.vnet.net
- Subject: [mg50990] Re: Simulating random fields in Mathematica
- From: Antti Penttilä <Antti.I.Penttila at helsinki.fi>
- Date: Fri, 1 Oct 2004 04:48:04 -0400 (EDT)
- Organization: University of Helsinki
- References: <cjgioq$q6v$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Mark Coleman wrote:
> Greetings,
>
> I'm searching for Mathematica code to simulate large (Gaussian) random fields,
> for work in spatial analysis. Can anyone point me in the right
> direction?
>
> Thanks,
>
> -Mark
>
I have a notebook for simulation of Gaussian random field as a Fourier-series representation. I have been told that this series represent Gaussian random field, but I don't have skills to theoretically check it myself. At least it seems to be almost normal.
Antti Penttilä
(************** Content-type: application/mathematica **************
CreatedBy='Mathematica 5.0'
Mathematica-Compatible Notebook
This notebook can be used with any Mathematica-compatible
application, such as Mathematica, MathReader or Publicon. The data
for the notebook starts with the line containing stars above.
To get the notebook into a Mathematica-compatible application, do
one of the following:
* Save the data starting with the line of stars above into a file
with a name ending in .nb, then open the file inside the
application;
* Copy the data starting with the line of stars above to the
clipboard, then use the Paste menu command inside the application.
Data for notebooks contains only printable 7-bit ASCII and can be
sent directly in email or through ftp in text mode. Newlines can be
CR, LF or CRLF (Unix, Macintosh or MS-DOS style).
NOTE: If you modify the data for this notebook not in a Mathematica-
compatible application, you must delete the line below containing
the word CacheID, otherwise Mathematica-compatible applications may
try to use invalid cache data.
For more information on notebooks and Mathematica-compatible
applications, contact Wolfram Research:
web: http://www.wolfram.com
email: info at wolfram.com
phone: +1-217-398-0700 (U.S.)
Notebook reader applications are available free of charge from
Wolfram Research.
*******************************************************************)
(*NotebookFileLineBreakTest
NotebookFileLineBreakTest*)
(*NotebookOptionsPosition[ 4509, 123]*)
(*NotebookOutlinePosition[ 5153, 145]*)
(* CellTagsIndexPosition[ 5109, 141]*)
(*WindowFrame->Normal*)
Notebook[{
Cell[CellGroupData[{
Cell["Gaussian surface as Fourier-series", "Title"],
Cell[BoxData[{
\(\(nn\ = \ 15;\)\ (*\
Length\ of\ the\ series\ *) \), "\[IndentingNewLine]",
\(\(\[Sigma]\ = \ 1;\)\ (*\
Deviation\ of\ the\ field\ *) \), "\[IndentingNewLine]",
\(\(\(covl\ = \ 2\)\(;\)\(\ \)\( (*\
The\ covariance\ length\ of\ the\ field\ \
*) \)\(\[IndentingNewLine]\)\(L\ = \ 20\)\(;\)\(\ \)\( (*\
The\ periodicity\ of\ the\ field\ will\ be\ [\(-L\), \
L]\ *) \)\)\)}], "Input"],
Cell[BoxData[
\(\(\( (*\
Misc . \ functions\ *) \)\(\[IndentingNewLine]\)\(sig2[k_,
n_]\ := \ \[Sigma]\^2\ \[Pi]\/2\ \((covl\/L)\)\^2\ Exp[\ \(-\(1\/2\)\
\)\ \((\(\[Pi]\ covl\)\/L)\)\^2\ \((k\^2 + n\^2)\)]\n
\[Phi][k_, n_]\ := \ Random[\ Real, \ {0, \ 2\ \[Pi]}]\n
A[k_, n_]\ := \ \@\(\(-2\)\ sig2[k, n]\ Log[\ Random[]]\)\)\)\)], "Input"],
Cell[CellGroupData[{
Cell[BoxData[
\(\(\( (*\
The\ series\ representation\ of\ the\ Gaussian\ field\ *) \)\(\
\[IndentingNewLine]\)\(z[x_, y_]\ := \
Sum[\ At[\([k + nn + 1,
n + nn +
1]\)]\ Cos[\ \(k\ \[Pi]\ x\)\/L\ + \ \(n\ \[Pi]\ y\)\/L\ + \
\ \[Phi]t[\([k + nn + 1, n + nn + 1]\)]], \ {k, \ \(-nn\), \
nn}, \ {n, \ \(-nn\), \ nn}]\)\)\)], "Input"],
Cell[BoxData[
RowBox[{\(General::"spell1"\), \(\(:\)\(\ \)\), "\<\"Possible spelling \
error: new symbol name \\\"\\!\\(\[Phi]t\\)\\\" is similar to existing symbol \
\\\"\\!\\(\[Phi]\\)\\\". \\!\\(\\*ButtonBox[\\\"More\[Ellipsis]\\\", \
ButtonStyle->\\\"RefGuideLinkText\\\", ButtonFrame->None, \
ButtonData:>\\\"General::spell1\\\"]\\)\"\>"}]], "Message"]
}, Open ]],
Cell[BoxData[
\(\(\( (*\ For\ a\ new\ realization\ of\ the\ field, \
calculate\ new\ coefficients\ for\ the\ series\ *) \)\(\
\[IndentingNewLine]\)\(\(\[Phi]t\ = \
Table[\ \[Phi][k, n], \ {k, \ \(-nn\), \ nn}, \ {n, \ \(-nn\), \
nn}];\)\[IndentingNewLine]
\(At\ = \
Table[\ A[k, n], \ {k, \ \(-nn\), \ nn}, \ {n, \ \(-nn\), \
nn}];\)\)\)\)], "Input"],
Cell[BoxData[
\(\(Plot3D[\ z[x, y], \ {x, \ \(-L\), \ L}, \ {y, \ \(-L\), \ L}, \
BoxRatios \[Rule] Automatic, \ PlotPoints \[Rule] 50];\)\)], "Input"],
Cell[BoxData[
\(\(\( (*\
Random\ sample\ *) \)\(\[IndentingNewLine]\)\(samp\ = \
Table[\ z[Random[Real, \ {\(-1\), 1}]\ L,
Random[Real, \ {\(-1\), 1}]\ L], \ {1000}];\)\)\)], "Input"],
Cell[BoxData[
\(\(\( (*\
Sample\ properties\ *) \)\(\[IndentingNewLine]\)\(Histogram[\ samp, \
HistogramScale \[Rule] 1];\)\)\)], "Input"],
Cell[BoxData[
\({Mean[\ samp], \ StandardDeviation[\ samp]}\)], "Input"]
}, Open ]]
},
FrontEndVersion->"5.0 for Microsoft Windows",
ScreenRectangle->{{0, 1280}, {0, 951}},
WindowSize->{1065, 923},
WindowMargins->{{0, Automatic}, {Automatic, 0}}
]
(*******************************************************************
End of Mathematica Notebook file.
*******************************************************************)