from IPython.core.display import HTML
css_file = 'pynoddy.css'
HTML(open(css_file, "r").read())

%matplotlib inline

# here the usual imports. If any of the imports fails, make sure that pynoddy is installed
# properly, ideally with 'python setup.py develop' or 'python setup.py install'
import sys, os
import matplotlib.pyplot as plt
import numpy as np
# adjust some settings for matplotlib
from matplotlib import rcParams
# print rcParams
rcParams['font.size'] = 15
# determine path of repository to set paths corretly below
repo_path = os.path.realpath('../..')
import pynoddy.history
import pynoddy.experiment
rcParams.update({'font.size': 20})

"""
Explanation: Topology analysis for the Gippsland Basin
In this notebook, we are testing the combination of the pynoddy Experiment class with the topology analysis and apply it to the Gippsland Basin model. This notebook will evaluate three aspects:

Application of the Experiment class to set-up a topology experiment
Running the topology analysis for a range of models
Evaluating and visualising results

Further steps might include the parallelisation of the random model generation, let's see where we get...
End of explanation
"""


reload(pynoddy.experiment)
# the model itself is now part of the repository, in the examples directory:
history_file = os.path.join(repo_path, "examples/GBasin_Ve1_V4_b.his")
gipps_topo_ex = pynoddy.experiment.Experiment(history = history_file) 

"""
Explanation: Creating an experiment object
First, we start with generating a pynoddy experiment object. The experiment class inherits all the methods from the base pynoddy.history class and we can directly import the Gippsland Basin model that we want to analyse into the object:
End of explanation
"""


gipps_topo_ex.info()

"""
Explanation: Some basic information about the model can be obtained with:
End of explanation
"""


gipps_topo_ex.plot_section('y')

"""
Explanation: We can have a quick look at the model in a section view (note that Noddy is now executed in the background when required - and the output automatically generated in the required resolution):
End of explanation
"""


# gipps_topo_ex.determine_model_stratigraphy()
gipps_topo_ex.plot_section('x', ve = 5, position = 'centre',
                cmap = 'YlOrRd',
                title = '',
                colorbar = False)
gipps_topo_ex.plot_section('y', position = 100, ve = 5.,
                cmap = 'YlOrRd',
                title = '',
                colorbar_orientation = 'horizontal')

"""
Explanation: The base plot is not very useful - but we can create a section plot with a define vertical exaggeration (keyword ve) and plot the colorbar in horizontal orientation:
End of explanation
"""


reload(pynoddy.experiment)
# the model itself is now part of the repository, in the examples directory:
history_file = os.path.join(repo_path, "examples/GBasin_Ve1_V4_b.his")
gipps_topo_ex = pynoddy.experiment.Experiment(history = history_file) 
gipps_topo_ex.load_parameter_file(os.path.join(repo_path, "examples/gipps_params.csv"))

"""
Explanation: Note: The names of the model stratigraphy (colorbar labels) are unfortunately not defined correctly in the input file - we need to fix that, then we should get useful labels, as well!
Loading parameters ranges from file
We now need to define the parameter ranges. This step can either be done through explicit definition in the notebook (see the previous notebook on the Experiment class), or a list of parameters and defined ranges plus statistics can be read in from a csv file. This enables the convenient parameter definition in a spreadsheed (for example through Excel).
In order to be read in correctly, the header should contain the labels:

'event' : event id
'parameter' : Noddy parameter ('Dip', 'Dip Direction', etc.)
'min' : minimum value
'max' : maximum value
'initial' : initial value

In addition, it is possible to define PDF type and parameters. For now, the following settings are supported:

'type' = 'normal' 
'stdev' : standard deviation
'mean' : mean value (default: 'initial' value)

We can read in the parameters simply with:
End of explanation
"""


gipps_topo_ex.freeze()

"""
Explanation: Generating random perturbations of the model
Before generating random prerturbations, we should now store the base version so that we can always revert to it at a later stage:
End of explanation
"""


gipps_topo_ex.set_random_seed(12345)

"""
Explanation: For a reproducible experiment, we can also set the random seed:
End of explanation
"""


gipps_topo_ex.random_perturbation()

"""
Explanation: And now, let's perturb the model:
End of explanation
"""


fig = plt.figure(figsize = (12,6))
ax1 = fig.add_subplot(211)
ax2 = fig.add_subplot(212)
gipps_topo_ex.plot_section(ax = ax1, direction = 'x', model_type = "base", 
                           colorbar = False, title = "", ve = 5.)
gipps_topo_ex.plot_section(ax = ax2, direction = 'x', colorbar = False, 
                           title = "", ve = 5.)

b1 = gipps_topo_ex.get_section('x', resolution = 50, model_type = 'base')
# b1.plot_section(direction = 'x', colorbar = False, title = "", ve = 5.)
b2 = gipps_topo_ex.get_section('x', resolution = 50, model_type = 'current')
diff = b1 - b2

print b1

b1.shape

"""
Explanation: Let's see what happened: we can compare the new model to the base model as we stored it before:
End of explanation
"""


gipps_topo_ex.random_perturbation()

fig = plt.figure(figsize = (12,6))
ax1 = fig.add_subplot(311)
ax2 = fig.add_subplot(312)
ax3 = fig.add_subplot(313)
gipps_topo_ex.plot_section(ax = ax1, direction = 'x', model_type = "base", 
                           colorbar = False, title = "", ve = 5.)
gipps_topo_ex.plot_section(ax = ax2, direction = 'x', colorbar = False, 
                           title = "", ve = 5.)
# plot difference


fig = plt.figure(figsize = (12,6))
ax1 = fig.add_subplot(211)
ax2 = fig.add_subplot(212)
gipps_topo_ex.plot_section(ax = ax1, direction = 'x', model_type = "base", 
                           colorbar = False, title = "", ve = 5.)
gipps_topo_ex.plot_section(ax = ax2, direction = 'x', colorbar = False, 
                           title = "", ve = 5.)

"""
Explanation: ...and another perturbation:
End of explanation
"""


