GU.calculate_entropy()

"""
Explanation: Information Entropy
As a first step, let's have a look at probabilities and then calculate information entropy for the entire block model!
End of explanation
"""


pickle.dump(GU.entropy, open("GB_entropy_high.pkl", "w"))
GU.export_to_vtk(GU.entropy, "entropy_high")

"""
Explanation: Save results and create VTK file for visualisations
End of explanation
"""


modes, bins = scipy.stats.mode(GU.all_blocks)
bins = bins[0,:,:,:]
strat_var_1 = modes[0,:,:,:]

"""
Explanation: Stratigraphic Variability
The modes can simply be calculated using the scipy stats mode method:
End of explanation
"""


s = bins.shape
strat_var_new = np.ndarray((s))
for i in range(s[0]):
    for j in range(s[1]):
        for k in range(s[2]):
            strat_var_new[i,j,k] = len(np.unique(GU.all_blocks[:,i,j,k]))


"""
Explanation: Actually, the correct way is to determine the number of unique values for stratigraphic variability! So: determine the number of unique elements (i.e. possible outputs) in each cell and combine with probability of mode, with bins calculated above:
End of explanation
"""


strat_var_2 = 1. - bins / GU.n
strat_combined = strat_var_new + strat_var_2

"""
Explanation: Thankfully, the second return value from the mode function is actually the bin count of the modal value, so exactly what we need for the second part of the stratigraphic variability function:
End of explanation
"""


pickle.dump(strat_var_1, open("strat_num_high.pkl", "w"))
pickle.dump(strat_var_2, open("strat_differs_high.pkl", "w"))
pickle.dump(strat_combined, open("strat_combined_high.pkl", "w"))

GU.export_to_vtk(strat_var_new, "strat_num_high")
GU.export_to_vtk(strat_var_2, "strat_bins_high")
GU.export_to_vtk(strat_combined, "strat_combined_high")

# Extract slice for information theory analysis
GU_slice = GU.all_blocks[:,0,:,:]
pickle.dump(GU_slice, open("Gippsland_x_slice.pkl", "w"))
GU_y_slice = GU.all_blocks[:,:,0,:]
pickle.dump(GU_y_slice, open("Gippsland_y_slice.pkl", "w"))


pwd

print a
print np.sort(a)
print np.argsort(a)
print "----"
print np.unique(a)
print np.argmax(np.unique(np.sort(a)))

a_sort = np.sort(a)
print a_sort

np.where(a_sort[:-1] != a_sort[1:])
n = float(len(a))

a = np.random.randint(0,3,size=(300,4))
b = [np.sum(a == id_a, axis=0) / 300. for id_a in np.unique(a)]

b

b

2**2 + 2**1 + 2**(-2) + 2**(-3)

def own_sqrt(x):
    s = 1.
    for k in range(4):
        s = 0.5 * (s + x/s)
        print s
    return s

own_sqrt(2.)

pynoddy?

cd ~/git/pynoddy/examples/


NH = pynoddy.NoddyHistory("GBasi123.his")

run gipps_uncertainty

for event in NH.events:
    print NH.events[event]

NH.events[2].event_type

type(NH.events)

NH.events

np.mod(101,10)

a = (1,2,3)
print np.random.choice(a, size = 3, replace=False)


"""
Explanation: Save results and create VTK visualisations:
End of explanation
"""


