You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
744 lines
135 KiB
Plaintext
744 lines
135 KiB
Plaintext
6 years ago
|
{
|
||
|
|
"cells": [
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 1,
|
||
|
|
"metadata": {
|
||
|
|
"ExecuteTime": {
|
||
|
|
"end_time": "2018-12-05T00:54:50.235522Z",
|
||
|
|
"start_time": "2018-12-05T00:54:42.731587Z"
|
||
|
|
}
|
||
|
|
},
|
||
|
|
"outputs": [],
|
||
|
|
"source": [
|
||
|
|
"import pandas as pd\n",
|
||
|
|
"import os\n",
|
||
|
|
"import numpy.ma as ma\n",
|
||
|
|
"\n",
|
||
|
|
"import numpy\n",
|
||
|
|
"from pyearth import Earth\n",
|
||
|
|
"from matplotlib import pyplot\n",
|
||
|
|
"\n",
|
||
|
|
"np.random.seed(2017)"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 2,
|
||
|
|
"metadata": {
|
||
|
|
"ExecuteTime": {
|
||
|
|
"end_time": "2018-12-05T00:54:54.936556Z",
|
||
|
|
"start_time": "2018-12-05T00:54:50.271465Z"
|
||
|
|
}
|
||
|
|
},
|
||
|
|
"outputs": [
|
||
|
|
{
|
||
|
|
"name": "stdout",
|
||
|
|
"output_type": "stream",
|
||
|
|
"text": [
|
||
|
|
"Importing profiles.csv\n"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"name": "stderr",
|
||
|
|
"output_type": "stream",
|
||
|
|
"text": [
|
||
|
|
"C:\\Users\\z5189959\\Desktop\\nsw-2016-storm-impact\\.venv\\lib\\site-packages\\numpy\\lib\\arraysetops.py:522: FutureWarning: elementwise comparison failed; returning scalar instead, but in the future will perform elementwise comparison\n",
|
||
|
|
" mask |= (ar1 == a)\n"
|
||
|
|
]
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"source": [
|
||
|
|
"def df_from_csv(csv, index_col, data_folder='../data/interim'):\n",
|
||
|
|
" print('Importing {}'.format(csv))\n",
|
||
|
|
" return pd.read_csv(os.path.join(data_folder,csv), index_col=index_col)\n",
|
||
|
|
"\n",
|
||
|
|
"df_profiles = df_from_csv('profiles.csv', index_col=[0, 1, 2])"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"metadata": {
|
||
|
|
"ExecuteTime": {
|
||
|
|
"end_time": "2018-12-04T23:49:04.770025Z",
|
||
|
|
"start_time": "2018-12-04T23:49:04.265699Z"
|
||
|
|
}
|
||
|
|
},
|
||
|
|
"source": [
|
||
|
|
"## Try using pyearth"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 73,
|
||
|
|
"metadata": {
|
||
|
|
"ExecuteTime": {
|
||
|
|
"end_time": "2018-12-05T01:54:00.320555Z",
|
||
|
|
"start_time": "2018-12-05T01:53:58.905803Z"
|
||
|
|
},
|
||
|
|
"code_folding": [
|
||
|
|
5,
|
||
|
|
20,
|
||
|
|
31,
|
||
|
|
40
|
||
|
|
],
|
||
|
|
"scrolled": false
|
||
|
|
},
|
||
|
|
"outputs": [
|
||
|
|
{
|
||
|
|
"data": {
|
||
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEWCAYAAABrDZDcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzs3Xl4ZGWZ8P/vXXtlreyd9J5uOr3TYDeNg4CAwoALig6LDuIy4juKzshwiSC+Oj+HURzGV0bHUUYRdQRkUVBEhh1UaKQbeqW7odNL9j2VVCW1n+f3Ry2k0kkn3Z2kstyf68qV1HlOznlOVXLu8+xijEEppdTcZct1BpRSSuWWBgKllJrjNBAopdQcp4FAKaXmOA0ESik1x2kgUEqpOU4DgZrzROQsEXlTRIIi8gER+YOIXJNK+7iI/CnXeVRqMmkgUMdFRL4uIv+T63xMsP8P+L4xpsAY87Ax5mJjzM8m84Qi8k4RaRpjn7tFxIjIGUO2LReRowb/pPaNi0jNsO1fF5FYKsj5ReRFEXn7sHxYqfSAiOwXkU+McHwRkYMi8voIaW4RuUtE+kWkTUSuH5Z+gYjsE5FBEXlWRBYPSbs8ladBEXluhGOfLyKvpo59UESuPdZ7pk6MBoI5SEQcc+Xc4zzfYmDPZOflBPUA/3KsHUQkH/gQ0Ad8dIRdfmWMKQDKgWeBB4alt6TSi4AvAv8tInXD9jkHqARqRWTTsLSvA6eQfB/PA74kIn+dyls58Gvgq0ApsBX41bDr+y7wrRGuywn8BvgRUAxcAXxHRE4d8Y1QJ0wDwSwhIodF5CYReV1EekXkpyLiSaW9U0SaRORGEWkDfpra/l4R2T7kSXH9kOPdKCLNQ54SL0j9c98MXJF6gtyR2rdGRH4rIj0ickBEPj3kOF8XkQdF5H9EpB/4eGrbA6ltARHZJSIrUvnvEJFGEblwgq/106m89aTyWpPaXg/UAr9LXZNbRJ4Tkb8b5dwrReTJ1HH2i8jlx8jnJ0Rkb+oaD4rIZ1Lb84E/ADWpcwaHP8kP8TNgvYicO9p5SAYBP8mSzTWj7WSMiQO/BOaLSMUI6cYY8xjJm/P6YcnXAI8Aj41wjo8B3zDG9Bpj9gL/DXw8lXYZsMcY84AxJkwyaJwqIitT53zKGHM/0DJClktJBqdfpPL2CrAXWD3aNaoTo4FgdvkocBGwDFgB3DIkbR7Jf6zFwLUicjpwF/AZoIzkU9dvUzfCOuA6YJMxpjB1zMPGmMeBfyX1hGmMST+Z3Qs0ATXAh4F/FZELhpz7UuBBwEfyRgTwPuAXQAnwGvC/JP8e55O8of1oAq/1fOCbwOVANXAEuA/AGLMMaADel7qmyGgnTN3AnwTuIfl0fBXwAxFZM8qvdADvJXkz+wTw/0TkdGPMAHAxqSfx1NdIN0KAQZLv+a3HeC+uIfkZ3AesTH22I+XfRfKm3Q30jpBuE5H3kyw5HBiyPY/k5/rL1NeVqWMhIiUkP/cdQw61A0i/J2uGpqWuvX5I+qiMMe2p6/qEiNglWaW1GNA2mwmmgWB2+b4xptEY00PyxnHVkDQL+JoxJmKMCQGfBn5kjHnZGJNI1YlHgDOBBOAGVouI0xhz2BhTP9IJRWQh8A7gRmNM2BizHfgxcPWQ3V5K1b1bqXMD/NEY87+pp9QHgArgW8aYGMkb2hIR8U3QtX4UuMsY82rqRn8T8HYRWXKM44/kvSQD4k+NMXFjzKvAQyRvkkcxxvzeGFOfepp9HngCOPs4zwnJoLhIRC4eniAii0hWx9yTunE+zdFP7JeLiB9If+4fTr3vaTVD0n8DXG+MeW1I+mUk/zaeAB4FHMB7UmkFqe99Q/bvAwqHpA9NG54+lnuB/5s6/x+BrxhjGsf5u2qcNBDMLkP/QY6QfFJL60wVzdMWA/+Uqhbyp24EC4EaY8wB4B9JFuM7ROS+Y1Rd1AA9xpjAsHPPHyVfae1Dfg4BXcaYxJDX8NZNZiTHc601qX0AMMYEST4VD83jeCwGNg97zz5KsgRyFBG5WES2pKqR/MAlJJ+2j0sqeH0j9SXDkq8G9qYCMCSf2D+Sql9Pu98Y4wOqgN3A24YdoyWVXgT8B3D+sPRrUseIp/Lya94KNsHU96Ih+xcBgSHpQ9OGp48qVX30K5KlGBfJUsSXROQ9x/xFddw0EMwuC4f8vIjsetfhPU0agVuNMb4hX3nGmHsBjDH3GGPeQfLmZ4DbRjlOC1AqIkOf8BYBzcc490Q4nmttIXkdQKaKp4zsPI5HI/D8sPeswBjz98N3FBE3ydLC7UBV6kb7GG/dyI/3PfkpyQbTDw7b/jGSDbhtqTaR75AMNkeVHowxXSSrAr8uItUjpEeAG4F1IvKB1HUsIBkY/nbIOT4MXCIi5caYXqAVGNqAeypvNb7vGZqWeu+XMb7G+bXA/lTJ0TLG7Ad+P9K1qZOjgWB2+ZyILBCRUpKNur86xr7/DfwfEdksSfki8h4RKRSROkl223MDYZJP6Omn9XaS1TY2gFQx/UXgmyLikWSD86d4qy1gshzPtd5Dsp55Q+qa/hV42Rhz+DjP+SiwQkSuFhFn6muTiKwaYV8Xyeq1TiCeqtYZ2gDeDpSJSPF4Tpyqyvk6yRs1AKk682XAGcCG1Ndaktc7YqOxMWYfyfaYL42SHgX+nWR1DCRLHG8AdUPOsYJkm1C6Ou7nwC0iUpJ6iv80cHcq7TfAWhH5UKpB//8CO1P5IFX37yFZ3WRL/Q2lSzOvAaek/hZFRJaRrJ4b2h6hJoAGgtnlHpL1uAdTX6N2OzTGbCX5D/t9kg2HB3irp4ebZHe+LqCNZMPozam0dNfDbhF5NfXzVcASkk/evyFZP//kRFzQMRzPtT5NsvviQySfXpcBVx7vCVPVXxemfreF5HtzG8n3a6R9vwDcT/L9/Qjw2yHp+0jWfx9MVTONVvU21L2p/KddAzxijNlljGlLfwF3AO9NBcmR/BvJRvTKUdLvItkm8b7UOX4w9Pipc/yQt4LN10g2AB8Bngf+LdWxAGNMJ8leTbem3ofNZL/3V5N80Pgvku0nIZIPKaTapT5JsrqqP3Xsh4CfjPE+qeMkRhemmRVE5DDwd8aYp3Kdl8k2l65VqamgJQKllJrjNBAopdQcp1VDSik1x2mJQCml5ricTT52PMrLy82SJUtynQ2llJpRtm3b1mWMOWpeqeFmRCBYsmQJW7duzXU2lFJqRhGRI2PvpVVDSik152kgUEqpOU4DgVJKzXEaCJRSao7TQKCUUnPcjOg1NNMZY4hEIsRiMZxOJ263G5Hh08orpVRuaCCYZJZl0dTUREdHBwBut5vKykoqKys1GCilpgUNBJPIsizq6+vZvn07lmVht9txOp10dnbicrnw+XwaDJRSOaeBYJIYY2hqamLHjh10d3cjIsRiMUQEp9NJIpFg3bp1VFVVaTBQSuWUNhZPkkgkQnt7O8YYXC4XALFYjEgkQjQapbW1lZ07d+L3+9GJ/5RSuaQlgkkSi8Ww2+04HA7sdjuhUIh4PA5AIpFgcHCQw4cPIyJaMlBK5ZQGgknicCTf2nQvoWAwiDEm005gWRaWZdHS0oLNZqOwsJD8/Pwc51opNRdpIJgExhj6+/szbQIABQUFGGMybQU2mw3LsgiFQhw4cACPx8OZZ56J3W7Pce6VUnONthFMgkgkQl9fH/Pnz2f16tWsW7cOp9NJa2srRUVFeDwejDF4vV4gGTgOHDjAG2+8gWVZOc69Umqu0UAwCdJP/DabDbfbTSgU4k9/+hNdXV1s27aNRCKB3W4nFothWRbGGKLRKLt37+bgwYMaDJRSU0oDwSRItwEYYwgGg9x7770kEgkAQqEQu3btYmBgAK/Xm6k6CoVCBAIBXn31VQ0GSqkppYFgErjdbnw+HwMDAwwODhKLxbLS0wPNjhw5QiKRIJFI4PV6McbQ19fH9u3bNRgopaaMBoJJICKUl5dTXl7Ok08+SSQSGXG/trY2GhoacLvdWJZFJBJhcHCQnp4etm/fTmNjo44xUEpNOg0Ek8AYQ2dnJ7/85S85cOBAVtr
|
||
|
|
"text/plain": [
|
||
|
|
"<Figure size 432x288 with 1 Axes>"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
"metadata": {
|
||
|
|
"needs_background": "light"
|
||
|
|
},
|
||
|
|
"output_type": "display_data"
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"name": "stdout",
|
||
|
|
"output_type": "stream",
|
||
|
|
"text": [
|
||
|
|
"[ { 'slope': -0.06617795750241562,\n",
|
||
|
|
" 'type': 'land',\n",
|
||
|
|
" 'x_end': 63.0,\n",
|
||
|
|
" 'x_start': 40.5,\n",
|
||
|
|
" 'z_end': 9.943990321404,\n",
|
||
|
|
" 'z_mean': 9.20564731619437,\n",
|
||
|
|
" 'z_start': 8.348161899068224},\n",
|
||
|
|
" { 'slope': -0.09424321140128608,\n",
|
||
|
|
" 'x_end': 69.0,\n",
|
||
|
|
" 'x_start': 63.0,\n",
|
||
|
|
" 'z_end': 10.532233612481347,\n",
|
||
|
|
" 'z_mean': 10.310033111599461,\n",
|
||
|
|
" 'z_start': 10.062895267627406},\n",
|
||
|
|
" { 'slope': 0.07278057230588836,\n",
|
||
|
|
" 'x_end': 96.0,\n",
|
||
|
|
" 'x_start': 69.0,\n",
|
||
|
|
" 'z_end': 9.071344271845113,\n",
|
||
|
|
" 'z_mean': 9.942154478828064,\n",
|
||
|
|
" 'z_start': 10.600758995514095},\n",
|
||
|
|
" { 'slope': 0.010222110457455531,\n",
|
||
|
|
" 'x_end': 105.0,\n",
|
||
|
|
" 'x_start': 96.0,\n",
|
||
|
|
" 'z_end': 8.940005362742282,\n",
|
||
|
|
" 'z_mean': 8.960960594412768,\n",
|
||
|
|
" 'z_start': 9.032076613769703},\n",
|
||
|
|
" { 'slope': -0.006783434031405405,\n",
|
||
|
|
" 'x_end': 171.0,\n",
|
||
|
|
" 'x_start': 105.0,\n",
|
||
|
|
" 'z_end': 9.377956662680647,\n",
|
||
|
|
" 'z_mean': 9.12742372282666,\n",
|
||
|
|
" 'z_start': 8.933411617025225},\n",
|
||
|
|
" { 'slope': 0.04487864391014401,\n",
|
||
|
|
" 'x_end': 177.0,\n",
|
||
|
|
" 'x_start': 171.0,\n",
|
||
|
|
" 'z_end': 9.159019876027982,\n",
|
||
|
|
" 'z_mean': 9.312931484934994,\n",
|
||
|
|
" 'z_start': 9.382296403384132},\n",
|
||
|
|
" { 'slope': 0.28077356578561746,\n",
|
||
|
|
" 'x_end': 198.0,\n",
|
||
|
|
" 'x_start': 177.0,\n",
|
||
|
|
" 'z_end': 3.9204468853264824,\n",
|
||
|
|
" 'z_mean': 6.663384763062234,\n",
|
||
|
|
" 'z_start': 9.046622396024778},\n",
|
||
|
|
" { 'slope': 0.08534897213040515,\n",
|
||
|
|
" 'x_end': 204.0,\n",
|
||
|
|
" 'x_start': 198.0,\n",
|
||
|
|
" 'z_end': 3.3242471595653202,\n",
|
||
|
|
" 'z_mean': 3.4936366209579663,\n",
|
||
|
|
" 'z_start': 3.754184455095742},\n",
|
||
|
|
" { 'slope': 0.02118839446027552,\n",
|
||
|
|
" 'type': 'berm',\n",
|
||
|
|
" 'x_end': 252.0,\n",
|
||
|
|
" 'x_start': 204.0,\n",
|
||
|
|
" 'z_end': 2.11417905440144,\n",
|
||
|
|
" 'z_mean': 2.7099232730105225,\n",
|
||
|
|
" 'z_start': 3.26349144192863},\n",
|
||
|
|
" { 'slope': 0.06273209072355844,\n",
|
||
|
|
" 'type': 'foreshore',\n",
|
||
|
|
" 'x_end': 285.5,\n",
|
||
|
|
" 'x_start': 252.0,\n",
|
||
|
|
" 'z_end': 0.0614306329519847,\n",
|
||
|
|
" 'z_mean': 1.0727411577535628,\n",
|
||
|
|
" 'z_start': 2.070753280242818}]\n"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"data": {
|
||
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEWCAYAAAB1xKBvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzs3XmcXFd16PvfqqGrep4ndWtyW5awJVvIssEXzGxjG7Cxr0NsCDFDcODCu0lu8gEywiOBS14SJyTcYDnEARJi42C4NsEBGxJsDEYgGU+yJEvqltTzXD3XvN4fdapc1V09SV1d3dXr+/n0p6r23lVn11HprLOHs4+oKsYYY8xiXPmugDHGmPXBAoYxxpglsYBhjDFmSSxgGGOMWRILGMYYY5bEAoYxxpglsYBhzBKJyGtE5ISITIrIO0XkP0TkDifvfSLyZL7raEwuWcAwOSEinxaRf8l3PVbYZ4AvqmqZqv5fVb1eVb+ayw2KyBtEpGuRMl8RERWRK9PSLhSRORdZOWWjIrJpVvqnRSTiBMOAiPxURK6aVY+4kz8hIsdF5P1ZPl9EpF1EXsyS5xORe0VkXET6ROR/zcp/s4gcE5FpEfkvEdmalvcup07TIvKjLJ/9DhF5wanfT0Xk4oX2mTk3FjDMvETEs1G2vcTtbQWO5Lou52gE+LOFCohIKfDfgTHgPVmKfENVy4A64L+Af5uV3+PkVwC/A/yDiOycVeZ1QANwgYhcMSvv08AOEvvxjcDHReQ6p251wLeAPwZqgEPAN2Z9v78BPp/le+0Avg58GKgCvgM8nM/fb6GygLHBiMhpEfl9EXlRREZF5J9ExO/kvUFEukTkEyLSB/yTk/52EXkm7czz0rTP+4SIdKeddb7ZOQj8AfCrzhnfs07ZTSLysIiMiMhJEflQ2ud8WkS+KSL/IiLjwPuctH9z0iZE5HkRucip/4CIdIrItSv8XT/k1G3EqesmJ/0UcAHwHec7+UTkRyLyG/Nse5eIPOZ8znERedcC9Xy/iBx1vmO7iPymk14K/Aewydnm5OyWQZqvApeKyOvn2w6JYBEg0VK6Y75CqholcQBuEZH6LPmqqo+QOIhfOiv7DuAh4JEs2/h14E9VdVRVjwL/ALzPybsFOKKq/6aqQRLB5TIR2eVs8weq+gDQk6XKbwV+rKpPOnX/c6AFWGhfmHNgAWNjeg+J/2RtwEXAH6XlNZE4w9sK3Cki+4B7gd8EaoEDJM7efM7Z5ceAK1S13PnM06r6PeBzOGesqnqZ89n3AV3AJuBW4HMi8ua0bd8EfJPEWeLXnbR3AP8MVAO/BL5P4nfbQuLAd2AFv+ubgP8NvAtoBs4A9wOoahtwFniH851C823QOdA/BvwribPt24G/F5FL5nnLAPB2Emfu7wf+WkT2qeoUcD3Omb3zl+2ACTBNYp9/doF9cQeJf4P7gV3Ov222+heROLgPA6NZ8l0iciOJlsjJtPQSEv+uX3f+bnM+CxGpJvHv/mzaRz0LJPfJJel5znc/lZa/EHH+Zr/evYT3mmWwgLExfVFVO1V1hMQB5va0vDjwKVUNqeoM8CHggKoeVNWY02cfAl4NxAAfcLGIeFX1tKqeyrZBEdkMvBb4hKoGVfUZ4MvAe9OKPeWMDcSdbUPizPH7zpnjvwH1wOdVNULiwLdNRKpW6Lu+B7hXVZ92AsLvA1eJyLYFPj+bt5MInP+kqlFVfRp4kMTBdA5V/a6qnnLO3B8HHgWuXuY2IRE8t4jI9bMzRGQLiW6gf1XVfuCHzG0BvEtEAkDy3/1WZ78nbUrL/zbwv1T1l2n5t5D4bTwK/DvgAd7m5JU5j2Np5ceA8rT89LzZ+Qt5DHi902osItG6LQJKlvBeswwWMDamzrTnZ0ic+SUNOl0CSVuB33W6owLOAWMzsElVTwK/TaL7YEBE7l+gy2QTMKKqE7O23TJPvZL6057PAEOqGkt7DS8fjLJZznfd5JQBQFUnSZxlp9dxKbYCr5q1z95DokUzh4hcLyI/c7qvAsANJM7el8UJcn/q/Mms7PcCR51ADYkWwLtFxJtW5gFVrQIagReAy2d9Ro+TXwH8LfCmWfl3OJ8RderyLV4OSpPOY0Va+QpgIi0/PW92/rxU9ZiznS8CvST23YskWrNmBVnA2Jg2pz3fQma/8OyZNZ3AZ1W1Ku2vRFXvA1DVf1XV15I4SCqJ/uNsn9MD1IhI+hnjFqB7gW2vhOV81x4S3wNIdS3VklnHpegEHp+1z8pU9SOzC4qIj0Tr4y+BRueA/AgvH/CXu0/+CagEbp6V/uskBqL7nDGbu0gcWOe0RlR1iEQX5KdFpDlLfgj4BLBHRN7pfI9WEgHk19K2cStwg4jUqeooiYP5ZWkfdRkvTyI4kp7n7Ps2ljjJQFW/qaq7VbUW+BSJf8dfLOW9ZuksYGxMHxWRVhGpIdF8/8YCZf8B+LCIvEoSSkXkbSJSLiI7ReRNzkEvSOKMP3n230+iu8gFoKqdwE+B/y0ifkkMnH+Ql8cqcmU53/VfgfeLyF7nO30OOKiqp5e5zX8HLhKR94qI1/m7QkRekaVsEYluvUEg6nQnpQ/k9wO1IlK5lA07XUifJnFAB0AS02PbgCuBvc7fbhLfN+vgt3PW/n3g4/Pkh4G/Av7ESXov8BKwM20bF5E4y092A34N+CMRqXYGsz8EfMXJ+zawW0T+uzMx4U+A55x6ICJuJ90DuJzfUKp1JCKXO2XqSXTNfSf5XrNyLGBsTP9Kop+53fmbdzqmqh4i8R/7iyQGQE/y8swWH4lpjkNAH4kB3j9w8pJTModF5Gnn+e3ANhJn8t8mMX7w2Ep8oQUs57v+kMS0zgdJnA23Abctd4NOt9u1znt7SOybPyexv7KV/Z/AAyT277uBh9Pyj5EYqG53urfm6/JLd59T/6Q7gIdU9XlV7Uv+AV8A3u4E02z+gsRkgIZ58u8lMWbyDmcbf5/++c427ubloPQpEgPZZ4DHgb9wJkigqoMkZnF91tkPryJz37+XxAnJl0iM78yQOJlJ+gKJGWDHnccPYVac2A2UNhYROQ38hqr+IN91ybWN9F2NWQ3WwjDGGLMkOQsYIrJZEpf3HxWRIyLyW056jXNB0wnnsXqe99/hlDkhzno9xhhj8idnXVLO7IpmVX3amRlzGHgnif7vEVX9vIh8EqhW1U/Mem9yaYD9JGaJHAYud2ZaGGOMyYOctTBUtde5YCk5sHeUxHz2m0gsY4Dz+M4sb38r8JiqjjhB4jHgulzV1RhjzOJWZXEu50rZVwIHScw174VEUJlnBkYLmRdcdTHPxVMicidwJ0Bpaenlu3btWrmKG2NMgTt8+PCQqs5ZMyybnAcMESkjMU3xt1V1XGT2BajZ35YlLWvfmareA9wDsH//fj106NC5VtUYYzYcETmzeKmEnM6Sci6seRD4uqp+y0nuT1496jwOZHlrF5lX6LaSfZVKY4wxqySXs6QE+EcS69fclZb1MC9fyJNcCnm27wPXOleEVpO4COr7uaqrMcaYxeWyhfEaEldnvkkS91J4RkRuIHFl8DUicgK4xnmNiOwXkS8DOCuL/imJtWB+AXzGSTPGGJMnBXWlt41hGGPM8ojIYVXdv5SydqW3McaYJbGAsc689NJLxGKxxQsaY8wKs5ukryPd3d3cddddNDU1ceONN/LKV76SJU5TNsaY82YBY42La5zOsU6ODR3j0fsfRVXp7e3lwIEDbN26lZtvvplXvCLbbRaMMWZl2aD3GhGOhTkxfIJjQ8c4OnSUo0NHOTZ0jGNDx5iOTNMYbuSmoZuyvnfnzp3cfPPNbN++fZVrbYxZ75Yz6G0BY5WNBcdSgSA9MJwaOUVM5x+bKI2Vsn98PztmduCaZ+jprW99K7fcckuuqm6MKUDLCRjWJZUDqkrvZC9HB4/OCQw9E+d2wfqUe4rHqx/n2bJn2T+xnwuCF8wps2PHDuLxOOPj44yOjhIOhwEoKyujurqa4uJiG/MwxpwzCxjnIRqP0j7anjUwjIfGV2w7Vf4qXlH3isRf/SvYVbeLqmAVP//hzzl
|
||
|
|
"text/plain": [
|
||
|
|
"<Figure size 432x288 with 1 Axes>"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
"metadata": {
|
||
|
|
"needs_background": "light"
|
||
|
|
},
|
||
|
|
"output_type": "display_data"
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"name": "stdout",
|
||
|
|
"output_type": "stream",
|
||
|
|
"text": [
|
||
|
|
"[ { 'slope': -0.006232531713848901,\n",
|
||
|
|
" 'type': 'land',\n",
|
||
|
|
" 'x_end': 87.0,\n",
|
||
|
|
" 'x_start': 54.0,\n",
|
||
|
|
" 'z_end': 18.907362973310256,\n",
|
||
|
|
" 'z_mean': 18.226859544431967,\n",
|
||
|
|
" 'z_start': 18.244707421270647},\n",
|
||
|
|
" { 'slope': 0.10048504940215892,\n",
|
||
|
|
" 'x_end': 90.0,\n",
|
||
|
|
" 'x_start': 87.0,\n",
|
||
|
|
" 'z_end': 18.81114066206634,\n",
|
||
|
|
" 'z_mean': 18.95124513461371,\n",
|
||
|
|
" 'z_start': 19.011103385573893},\n",
|
||
|
|
" { 'slope': 0.7867562476656698,\n",
|
||
|
|
" 'x_end': 105.0,\n",
|
||
|
|
" 'x_start': 90.0,\n",
|
||
|
|
" 'z_end': 8.264119331473214,\n",
|
||
|
|
" 'z_mean': 13.443236505259332,\n",
|
||
|
|
" 'z_start': 18.482680684877423},\n",
|
||
|
|
" { 'slope': 0.138231538826561,\n",
|
||
|
|
" 'x_end': 108.0,\n",
|
||
|
|
" 'x_start': 105.0,\n",
|
||
|
|
" 'z_end': 7.625424767531349,\n",
|
||
|
|
" 'z_mean': 7.719382219664138,\n",
|
||
|
|
" 'z_start': 7.899377829775406},\n",
|
||
|
|
" { 'slope': -0.24669278691155727,\n",
|
||
|
|
" 'x_end': 111.0,\n",
|
||
|
|
" 'x_start': 108.0,\n",
|
||
|
|
" 'z_end': 8.25183247996525,\n",
|
||
|
|
" 'z_mean': 8.037437899312298,\n",
|
||
|
|
" 'z_start': 7.767262334494268},\n",
|
||
|
|
" { 'slope': -0.022224358740496428,\n",
|
||
|
|
" 'x_end': 138.0,\n",
|
||
|
|
" 'x_start': 111.0,\n",
|
||
|
|
" 'z_end': 8.805851879663944,\n",
|
||
|
|
" 'z_mean': 8.69776858676918,\n",
|
||
|
|
" 'z_start': 8.321441165942321},\n",
|
||
|
|
" { 'slope': 0.011690412029154119,\n",
|
||
|
|
" 'x_end': 174.0,\n",
|
||
|
|
" 'x_start': 138.0,\n",
|
||
|
|
" 'z_end': 8.399119455457473,\n",
|
||
|
|
" 'z_mean': 8.591198407135863,\n",
|
||
|
|
" 'z_start': 8.769268316237742},\n",
|
||
|
|
" { 'slope': 0.06596545786924302,\n",
|
||
|
|
" 'x_end': 183.0,\n",
|
||
|
|
" 'x_start': 174.0,\n",
|
||
|
|
" 'z_end': 7.84469125226732,\n",
|
||
|
|
" 'z_mean': 8.148131103000907,\n",
|
||
|
|
" 'z_start': 8.376950179198564},\n",
|
||
|
|
" { 'slope': 0.19014969366205653,\n",
|
||
|
|
" 'x_end': 213.0,\n",
|
||
|
|
" 'x_start': 183.0,\n",
|
||
|
|
" 'z_end': 2.9990666738937883,\n",
|
||
|
|
" 'z_mean': 5.087790707786215,\n",
|
||
|
|
" 'z_start': 7.748699573256347},\n",
|
||
|
|
" { 'slope': 0.012096582737314357,\n",
|
||
|
|
" 'type': 'berm',\n",
|
||
|
|
" 'x_end': 219.0,\n",
|
||
|
|
" 'x_start': 213.0,\n",
|
||
|
|
" 'z_end': 2.9163517442302123,\n",
|
||
|
|
" 'z_mean': 2.9347601110839676,\n",
|
||
|
|
" 'z_start': 2.975348176164513},\n",
|
||
|
|
" { 'slope': 0.03056014910314705,\n",
|
||
|
|
" 'x_end': 255.0,\n",
|
||
|
|
" 'x_start': 219.0,\n",
|
||
|
|
" 'z_end': 2.050738949067936,\n",
|
||
|
|
" 'z_mean': 2.653057687157337,\n",
|
||
|
|
" 'z_start': 2.927273834441084},\n",
|
||
|
|
" { 'slope': 0.06806395674070564,\n",
|
||
|
|
" 'type': 'foreshore',\n",
|
||
|
|
" 'x_end': 285.5,\n",
|
||
|
|
" 'x_start': 255.0,\n",
|
||
|
|
" 'z_end': 0.04691924522162938,\n",
|
||
|
|
" 'z_mean': 1.100297110086715,\n",
|
||
|
|
" 'z_start': 2.0068617905167647}]\n"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"data": {
|
||
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAEWCAYAAABliCz2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzs3XmcZGV96P/Pc2rfe63ee6ZnY2BYZBhBRCIiUdC4xKuJmngNmpDN5N7f/eW63V9cck2MiVlMNIkkgjFuMUaMKIKgDAiyyDIwDMzQMz0zve9d3V37qTrP74+qOnT13j1d09v3/Xr1i65znjrnqRr6W099n01prRFCCLH1GetdASGEEOeHBHwhhNgmJOALIcQ2IQFfCCG2CQn4QgixTUjAF0KIbUICvtjWlFLXKKU6lVJxpdRblVI/VEq9t3juN5RSD613HYVYKxLwtzml1JeVUp9aZtlPKKW+Wuk6nWd/Anxeax3UWn9Xa32T1vpfK3lDpdR1SqneZZS7Uil1l1IqppQaV0o9rpS6ecY1rOIH1bRS6kTp3Izna6VUolim9PPB4rkqpdRtSqnB4vNfVEp9aAWv4d1KqbPF639XKVUz41yNUuqO4rmzSql3zzjXpJT6nlKqv1i/nbOue2xWfXNKqTuXWy+xOAn44rxRSjk34P12AMcqXZeVUkpdDfwEeADYA9QCvwvcNKNYv9Y6CISB/wf4Z6XUBbMudVnxw6z08xfF438DBIELgQjwZuDUMut2APgi8B6gAUgC/zCjyBeAbPHcrwH/WHwOgAXcDfy3+a6ttT5QqisQArqB/1hOvcQyaK3lZ5P8AGeAjwDPAxPA7YB3xvnfAk4C48D3gObicUXhD3wYmASeBS4GbgFMCn+cceDOYvkPAX3ANHACeC1wY7GcWSz7TLFsc/Fe48V7/9aM+nwC+DbwVWAK+M3isf8oHpsGjgL7iq9rGOgBXrea9wC4Dugt1n8Q+Lcl3pdTFAJQqviaPMBh4DeL538DeGjGvfcD9xavcwL4lUXqeTPwQvE1dgG/XTweKN7PKt4zXqrPrOc/BHxhketfB/TOOjYMvGPGYw3sWeD5zwFvXeT6GvgdoLP4Pn8BUMVzfwZ8fUbZ3cX/N0LF15cF9s04/2/An8+6vrN4j52L1OHVxfcnsN5/e1vlZ90rID8r+McqBLvngDagBngY+FTx3PXAKHCwGLj+HniweO71wJNAFYXgfyHQVDz35dI1io8vKAbdUlDcCewu/v4J4Kuz6vQAhdadF3gZMAK8dkZ5E3grhW+TvuKxdLFOTuArwGng/wAuCsH59Crfg+uAHPCZ4nvgW+x9mXG9G2Y8Psw8Ab8YyHooBHJn8XqjwIEF6vnGYiBUxcCVBA7OqGfvIq/RD+SB1yxSxr5G8b19M4UPkctnlFks4P8LhW82NwN75zmvge8X/59pL/673lg891/Ah2aVjwNXAJcDqVnn/ohiY2LGseUE/NuAL6/3391W+pGUzubzea11j9Z6HPhT4F3F478G3Ka1fkprnaHQCr66mCM1KbS+9lNopb2gtR5Y4Pp5CoHxIqWUS2t9Rms971d9pVQb8CoKf/xprfURCoHkPTOKPaILuXFLa50qHvup1voerXWOQmu/nkIL0AS+CexUSlWt4j2AQtD7uNY6U7zfYu/LSvwScEZrfbvWOqe1fgr4T+Dt8xXWWv9Aa31KFzwA/Ai4dpn3qqYQxBf6NyppVkrFKHxjuAP4X1rrp2eVearYB1D6eX3x+B8AXwM+ADyvlDqplLpp1nP/XGsd01p3A/dT+ECHQipoclbZSQr/jy12btmUUn4K7+2XV/I8sTgJ+JtPz4zfz1JIqVD879nSCa11HBgDWrTWPwE+T+Fr+ZBS6lalVHi+i2utTwL/k0JLfFgp9U2lVPN8ZYv3HNdaT8+qU8sC9S0ZmvF7ChjVWudnPIZC4FjIQu8BwIjWOj2rjvO+L4tcfz47gKtmBk8KHyaN8xVWSt2klHq02NkaA94A1C3zXhMUPrialijXr7WuopDD/zsK32ZmO6i1rprxcw+A1jqltf4zrfUVFPoHvgX8x8zOVwppsZIkL/2bxIv3nClMIX212LmVeBuF1NkDK3yeWIQE/M2nbcbv7UB/8fd+CkEJAKVUgMIfch+A1vrvin/cByjkzP93seic5VK11l/XWr+qeD1NIUUyX9l+oEYpNbP11l6650LXXwMLvQfz3W/R92UFeoAHZgXPoNb6d2cXVEp5KLT+Pws0FIPyXRTSO/PVsYzWOgk8wgIdm/OUz1Dot7hEKfXWZb+il54/RSEvHwA6lvGUY8BlpQdKqV0UvhW+WPxxKqX2zih/GSvvGH8v8BWttSznu4Yk4G8+v6+Uai22xD4K/Hvx+NeBm5VSLysGnD8DHtNan1FKvVwpdZVSygUkKOTQSy3qIWBX6eJKqQuUUtcXr5Gm0OKeWXanUsoA0Fr3AD8DPq2U8iqlLgXeTyFVUEkLvQfzWfB9WeE9vw/sU0q9RynlKv68XCl14Txl3RQC4AiQK6ZKXjfj/BBQq5SKLHK/DwK/oZT630qpWgCl1GVKqW/OV1hrnQX+CvjYcl6MUuqPi/V3K6W8wP8AYhQ6o5fyNeBNSqlrix+gfwJ8R2s9rbVOAN8B/kQpFVBKXQO8hULHbeneXgrvD4Cn+Hhm3VqB1wAVHR67HUnA33y+TiEf3FX8+RSA1vrHwB9TaFkOUOgwfGfxOWHgnymkCs5SSGl8tnjuSxTy9TGl1Hcp/CH+OYUOyUEgSiGowkvD48aUUk8Vf38XhY7dfgp55I9rre9d01c817zvwXyWeF+WrZi2el3xuf0U3ptS5/B8Zf+QQppkAng3hdFBpfPHgW8AXcX3fU7KTGv9MwopmuuL5caBWyl8U1jIbUC7UupNM449M2tc+9+WbkFhhNNo8fX8IvDGYsprqffiGIURPF+jMDIoBPzejCK/R6HDfLj4On+3+JyS0qgogOO8lMYreQ+Fvp9lDRMVy6fkG9PmoZQ6Q2EEyX3rXZf1Iu+BEKsnLXwhhNgmJOALIcQ2ISkdIYTYJqSFL4QQ28R5XcxqKXV1dXrnzp3rXQ0hhNg0nnzyyVGtdf1yym6ogL9z506eeOKJ9a6GEEJsGkqps0uXKpCUjhBCbBMS8IUQYpuQgC+EENuEBHwhhNgmJOALIcQ2saFG6QghxExaazKZDKZp4nK58Hg8KKWWfqKYlwR8IcSGpLVmZGSEiYkJ8vk8pmlSU1NDS0sLhiHJidWQgC+E2JAymQwTExOYpkk8XlhN+eTJkwC0trZKS38V5GNSCLEhZbNZUqkUY2NjOJ1OfD4fPp+P8fFxMpnMeldvU5IWvhBiw9FaMzk5yeDgIKZpkk6n8fv9OJ1OnE4npmni9XqXvpAoIy18IcSGk8lkSCaTNDQ0AGBZFqOjo2itsSwLp1PaqqshAV8Ise4OHz7M5OSk/dg0TRwOB/X19bS2tqK1xjRNpqenSaVSTE1NIUu7r5x8TAohKq40vDKbzaK1RimF2+3G4/Fw6tQpvvGNb/Dtb3+b66+/nte//vU4nU7S6TRaa/x+P8FgEL/fT3NzM36/n8nJSSKRiKR1VkgCvhCiokrDK8fHx5mYmGBqaopQKERtbS1VVVXccccdQKFVf88999Db28uv/uqvkkwmGRsbI51OY5omO3bswOVykU6nyeVyZLNZCfgrJAFfCFFR6XSaoaEhpqenGR4exuPxMDQ0BEBnZ6c91LLkta99LZOTkzQ3N5PL5UgkEvT19WGapv28VCpFKBQiFArJ8MwVkIAvhKgYrTUDAwP09fUxOTmJaZqYpkkul6O3t3fO/he7du2ivb2doaEhDMPA7XbjcrmYmppiZGSESCQCQH19PclkkkwmI638FZCAL4SomHQ6zcDAAIlEgnw+TzabJR6P4/P5iMViZR21AD09PTgcDvL5PJlMhnw+j8PhwOv14vf7CQQCOJ1OXC4XqVRKhmeukAR8IUTFxONxkskkwWCQVCqFYRik02kcDgdnzpyZU940TSYnJ9Fac+bMGTv419bW4vF
|
||
|
|
"text/plain": [
|
||
|
|
"<Figure size 432x288 with 1 Axes>"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
"metadata": {
|
||
|
|
"needs_background": "light"
|
||
|
|
},
|
||
|
|
"output_type": "display_data"
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"name": "stdout",
|
||
|
|
"output_type": "stream",
|
||
|
|
"text": [
|
||
|
|
"[ { 'slope': 0.0008328773436106309,\n",
|
||
|
|
" 'type': 'berm',\n",
|
||
|
|
" 'x_end': 58.0,\n",
|
||
|
|
" 'x_start': 37.0,\n",
|
||
|
|
" 'z_end': 1.8572595562110041,\n",
|
||
|
|
" 'z_mean': 1.8943805115420846,\n",
|
||
|
|
" 'z_start': 1.910394386927169},\n",
|
||
|
|
" { 'slope': -0.021220543877635122,\n",
|
||
|
|
" 'type': 'berm',\n",
|
||
|
|
" 'x_end': 138.0,\n",
|
||
|
|
" 'x_start': 58.0,\n",
|
||
|
|
" 'z_end': 3.5977328867797485,\n",
|
||
|
|
" 'z_mean': 2.7247048500781106,\n",
|
||
|
|
" 'z_start': 1.8468598845719795},\n",
|
||
|
|
" { 'slope': -0.06381773295732189,\n",
|
||
|
|
" 'x_end': 142.0,\n",
|
||
|
|
" 'x_start': 138.0,\n",
|
||
|
|
" 'z_end': 3.8359998018824712,\n",
|
||
|
|
" 'z_mean': 3.741723519647647,\n",
|
||
|
|
" 'z_start': 3.6479493664613014},\n",
|
||
|
|
" { 'slope': -0.016157748520191067,\n",
|
||
|
|
" 'x_end': 157.0,\n",
|
||
|
|
" 'x_start': 142.0,\n",
|
||
|
|
" 'z_end': 4.169920934351783,\n",
|
||
|
|
" 'z_mean': 4.013083138991948,\n",
|
||
|
|
" 'z_start': 3.8843177342116757},\n",
|
||
|
|
" { 'slope': -0.1322171985012303,\n",
|
||
|
|
" 'x_end': 177.0,\n",
|
||
|
|
" 'x_start': 157.0,\n",
|
||
|
|
" 'z_end': 6.427354838550462,\n",
|
||
|
|
" 'z_mean': 5.4360212510513115,\n",
|
||
|
|
" 'z_start': 4.223480595415142},\n",
|
||
|
|
" { 'slope': 0.34782985777944603,\n",
|
||
|
|
" 'x_end': 193.0,\n",
|
||
|
|
" 'x_start': 177.0,\n",
|
||
|
|
" 'z_end': 1.3683244048034364,\n",
|
||
|
|
" 'z_mean': 3.876447007928749,\n",
|
||
|
|
" 'z_start': 6.295019014555152},\n",
|
||
|
|
" { 'slope': 0.06019958288565148,\n",
|
||
|
|
" 'type': 'foreshore',\n",
|
||
|
|
" 'x_end': 213.5,\n",
|
||
|
|
" 'x_start': 193.0,\n",
|
||
|
|
" 'z_end': 0.06618164722605435,\n",
|
||
|
|
" 'z_mean': 0.6722026822899186,\n",
|
||
|
|
" 'z_start': 1.1980470061124415}]\n"
|
||
|
|
]
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"source": [
|
||
|
|
"from scipy.signal import savgol_filter\n",
|
||
|
|
"import re\n",
|
||
|
|
"from scipy.stats import linregress\n",
|
||
|
|
"import warnings\n",
|
||
|
|
"warnings.simplefilter(action='ignore', category=FutureWarning)\n",
|
||
|
|
"\n",
|
||
|
|
"def get_breakpoints(model, min_distance=20):\n",
|
||
|
|
" # Get breakpoints\n",
|
||
|
|
" breakpoints = []\n",
|
||
|
|
" for line in model.summary().split('\\n'):\n",
|
||
|
|
" # Get unpruned lines\n",
|
||
|
|
" if 'No' in line and 'None' not in line:\n",
|
||
|
|
" # Get break points\n",
|
||
|
|
" m = re.search(\"h\\(x0-(\\d*\\.?\\d+)\\)\", line)\n",
|
||
|
|
" if m:\n",
|
||
|
|
" breakpoints.append(float(m.groups()[0]))\n",
|
||
|
|
" m = re.search(\"h\\((\\d*\\.?\\d+)-x0\\)\", line)\n",
|
||
|
|
" if m:\n",
|
||
|
|
" breakpoints.append(float(m.groups()[0]))\n",
|
||
|
|
" return sorted(list(set(breakpoints)))\n",
|
||
|
|
" \n",
|
||
|
|
"def get_segments(breakpoints, x_min, x_max):\n",
|
||
|
|
" segments = []\n",
|
||
|
|
" breakpoints = [x_min] + breakpoints + [x_max]\n",
|
||
|
|
"\n",
|
||
|
|
" for x1, x2 in zip(breakpoints, breakpoints[1:]):\n",
|
||
|
|
" segments.append({\n",
|
||
|
|
" 'x_start': x1,\n",
|
||
|
|
" 'x_end': x2\n",
|
||
|
|
" })\n",
|
||
|
|
" return segments \n",
|
||
|
|
"\n",
|
||
|
|
"def get_segment_slopes(segments, x, z):\n",
|
||
|
|
" for segment in segments:\n",
|
||
|
|
" mask = ma.masked_where((segment['x_start'] < x) & (x < segment['x_end']),x ).mask\n",
|
||
|
|
" segment['z_mean'] = np.mean(z[mask])\n",
|
||
|
|
" segment['z_start'] = np.mean(z[mask][0])\n",
|
||
|
|
" segment['z_end'] = np.mean(z[mask][-1])\n",
|
||
|
|
" segment['slope'] = -linregress(x[mask], z[mask]).slope\n",
|
||
|
|
" return segments\n",
|
||
|
|
" \n",
|
||
|
|
"def classify_segments(segments, x,z):\n",
|
||
|
|
" \n",
|
||
|
|
" # Most seaward slope must be foreshore\n",
|
||
|
|
" segments[-1]['type'] = 'foreshore'\n",
|
||
|
|
" \n",
|
||
|
|
" # Most landward slope must be land\n",
|
||
|
|
" segments[0]['type'] = 'land'\n",
|
||
|
|
" \n",
|
||
|
|
" # Segments with really high slopes must be structures\n",
|
||
|
|
" for seg in segments:\n",
|
||
|
|
" if seg['slope'] > 2.0:\n",
|
||
|
|
" seg['type'] = 'structure'\n",
|
||
|
|
" \n",
|
||
|
|
" # Segments with large change of slope and \n",
|
||
|
|
" # Segment with max slope should be dune face\n",
|
||
|
|
"# dune_face_idx = [n for n, seg in enumerate(segments) if seg['slope']==max(x['slope'] for x in segments)][0]\n",
|
||
|
|
"# segments[dune_face_idx]['type'] = 'dune_face'\n",
|
||
|
|
" \n",
|
||
|
|
" # Pick out berms \n",
|
||
|
|
" for seg in segments:\n",
|
||
|
|
" if (-0.03 < seg['slope'] < 0.03 # berms should be relatively flat\n",
|
||
|
|
" and 0 < seg['z_mean'] < 4 # berms should be located between 0-4 m AHD\n",
|
||
|
|
" ): # berms should be seaward of dune face\n",
|
||
|
|
" seg['type'] = 'berm'\n",
|
||
|
|
" \n",
|
||
|
|
"# slope = None\n",
|
||
|
|
"# for seg in reversed(segments):\n",
|
||
|
|
"# if slope is None:\n",
|
||
|
|
"# continue\n",
|
||
|
|
"# elif slope - 0.03 < seg['slope'] < slope + 0.03:\n",
|
||
|
|
"# seg['type'] = 'foreshore'\n",
|
||
|
|
"# else:\n",
|
||
|
|
"# break\n",
|
||
|
|
" \n",
|
||
|
|
" return segments\n",
|
||
|
|
"\n",
|
||
|
|
"def get_piecewise_linear_model(x,z):\n",
|
||
|
|
" #Fit an Earth model\n",
|
||
|
|
" model = Earth(penalty=3,thresh=0.0005)\n",
|
||
|
|
" model.fit(x,z)\n",
|
||
|
|
" return model\n",
|
||
|
|
"\n",
|
||
|
|
"def plot_profile_classification(site_id, profile_type):\n",
|
||
|
|
" df_profile = df_profiles.query(\"site_id == '{}' and profile_type == '{}'\".format(site_id, profile_type))\n",
|
||
|
|
" x = np.array(df_profile.index.get_level_values('x').tolist())\n",
|
||
|
|
" z = np.array(df_profile.z.tolist()) \n",
|
||
|
|
" \n",
|
||
|
|
" nan_mask = ma.masked_invalid(z).mask\n",
|
||
|
|
" x = x[~nan_mask]\n",
|
||
|
|
" z_unfiltered = z[~nan_mask]\n",
|
||
|
|
" z = savgol_filter(z_unfiltered, 51, 3)\n",
|
||
|
|
" \n",
|
||
|
|
" model = get_piecewise_linear_model(x,z)\n",
|
||
|
|
" breakpoints = get_breakpoints(model)\n",
|
||
|
|
" segments = get_segments(breakpoints, x_min=x.min(), x_max=x.max())\n",
|
||
|
|
" segments = get_segment_slopes(segments, x=x, z=z)\n",
|
||
|
|
"# segments = merge_similar_segments(segments)\n",
|
||
|
|
" segments = classify_segments(segments, x=x, z=z)\n",
|
||
|
|
" \n",
|
||
|
|
" pyplot.figure()\n",
|
||
|
|
" pyplot.plot(x,z_unfiltered, color='0.5',marker='.', alpha=.2, ms=10,linestyle=\"None\")\n",
|
||
|
|
"\n",
|
||
|
|
" # Plot different segments\n",
|
||
|
|
" foreshore_segments = [x for x in segments if x.get('type') == 'foreshore']\n",
|
||
|
|
" for seg in foreshore_segments:\n",
|
||
|
|
" pyplot.plot([seg['x_start'], seg['x_end']],\n",
|
||
|
|
" [seg['z_start'], seg['z_end']],\n",
|
||
|
|
" linewidth=4, \n",
|
||
|
|
" color='b')\n",
|
||
|
|
"\n",
|
||
|
|
" land_segments = [x for x in segments if x.get('type') == 'land']\n",
|
||
|
|
" for seg in land_segments:\n",
|
||
|
|
" pyplot.plot([seg['x_start'], seg['x_end']],\n",
|
||
|
|
" [seg['z_start'], seg['z_end']],\n",
|
||
|
|
" linewidth=4, \n",
|
||
|
|
" color='g')\n",
|
||
|
|
"\n",
|
||
|
|
" berm_segments = [x for x in segments if x.get('type') == 'berm']\n",
|
||
|
|
" for seg in berm_segments:\n",
|
||
|
|
" pyplot.plot([seg['x_start'], seg['x_end']],\n",
|
||
|
|
" [seg['z_start'], seg['z_end']],\n",
|
||
|
|
" linewidth=4, \n",
|
||
|
|
" color='y')\n",
|
||
|
|
"\n",
|
||
|
|
" dune_face_segments = [x for x in segments if x.get('type') == 'dune_face']\n",
|
||
|
|
" for seg in dune_face_segments:\n",
|
||
|
|
" pyplot.plot([seg['x_start'], seg['x_end']],\n",
|
||
|
|
" [seg['z_start'], seg['z_end']],\n",
|
||
|
|
" linewidth=4, \n",
|
||
|
|
" color='r')\n",
|
||
|
|
" \n",
|
||
|
|
" structure_segments = [x for x in segments if x.get('type') == 'structure']\n",
|
||
|
|
" for seg in structure_segments:\n",
|
||
|
|
" pyplot.plot([seg['x_start'], seg['x_end']],\n",
|
||
|
|
" [seg['z_start'], seg['z_end']],\n",
|
||
|
|
" linewidth=4, \n",
|
||
|
|
" color='m')\n",
|
||
|
|
" \n",
|
||
|
|
" unclassified_segments = [x for x in segments if x.get('type') is None]\n",
|
||
|
|
" for seg in unclassified_segments:\n",
|
||
|
|
" pyplot.plot([seg['x_start'], seg['x_end']],\n",
|
||
|
|
" [seg['z_start'], seg['z_end']],\n",
|
||
|
|
" linewidth=4, \n",
|
||
|
|
" color='0.4')\n",
|
||
|
|
"\n",
|
||
|
|
" pyplot.xlabel('x (m)')\n",
|
||
|
|
" pyplot.ylabel('z (m AHD)')\n",
|
||
|
|
" pyplot.title('{} profile at {}'.format(profile_type, site_id))\n",
|
||
|
|
" pyplot.show()\n",
|
||
|
|
"\n",
|
||
|
|
" import pprint\n",
|
||
|
|
" pp = pprint.PrettyPrinter(indent=4)\n",
|
||
|
|
" pp.pprint(segments)\n",
|
||
|
|
"\n",
|
||
|
|
"plot_profile_classification('NARRA0018', 'prestorm')\n",
|
||
|
|
"plot_profile_classification('NARRA0019', 'prestorm')\n",
|
||
|
|
"plot_profile_classification('CRESn0017', 'poststorm')"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"metadata": {
|
||
|
|
"heading_collapsed": true
|
||
|
|
},
|
||
|
|
"source": [
|
||
|
|
"## Try lmfit"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": null,
|
||
|
|
"metadata": {
|
||
|
|
"hidden": true
|
||
|
|
},
|
||
|
|
"outputs": [],
|
||
|
|
"source": []
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 26,
|
||
|
|
"metadata": {
|
||
|
|
"ExecuteTime": {
|
||
|
|
"end_time": "2018-12-05T01:13:41.773484Z",
|
||
|
|
"start_time": "2018-12-05T01:13:27.230746Z"
|
||
|
|
},
|
||
|
|
"code_folding": [
|
||
|
|
0
|
||
|
|
],
|
||
|
|
"hidden": true,
|
||
|
|
"scrolled": true
|
||
|
|
},
|
||
|
|
"outputs": [
|
||
|
|
{
|
||
|
|
"data": {
|
||
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD8CAYAAABn919SAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3Xd4XNd95//3uXd6AQYdIApBgCRYQIqkQJGy4shFSizbWTt5fonlEstZ29okdhLncVm3NG82VrL7s39pdlZx3GRbttcltmNv1ip2LFsSxQ6CnQTRiF6mYPrMPb8/BnPFLpIAUYbf1/PgIQWOMOcC5GfOfO/3nKO01gghhFj5jKUegBBCiIUhgS6EECVCAl0IIUqEBLoQQpQICXQhhCgREuhCCFEiJNCFEKJESKALIUSJkEAXQogS4VjMJ6uurtatra2L+ZRCCLHi7d+/f1JrXfNSj1vUQG9tbWXfvn2L+ZRCCLHiKaX6r+dxUnIRQogSIYEuhBAlQgJdCCFKhAS6EEKUCAl0IYQoEYva5bLSaa1Jp9Nks1mcTidutxul1FIPSwghAAn062ZZFkNDQ4yPjwPgdrupra2ltrZWQl0IsSxIoF8Hy7I4e/Yshw4dwrIsTNPE6XQyMTGBy+UiFApJqAshlpwE+kvQWjM0NMThw4eZmppCKUU2m0UphdPpJJ/Ps2XLFurq6iTUhRBLSm6KvoR0Os3Y2Bhaa1wuFwDZbJZ0Ok0mk2FkZITu7m7C4TBy4LYQYinJDP0lZLNZTNPE4XBgmibJZJJcLgdAPp8nkUjQ19eHUkpm6kKIJSWB/hIcjsK3qNjVMjs7i9barqNbloVlWQwPD2MYBsFgEL/fv8SjFkLcjiTQr0FrTTQatWvmAIFAwP7zbDaLYRhYlkUymeTMmTN4PB52796NaZpLNWwhxG1KAv0a0uk0kUiExsZGamtrSSQSzM7O4vF4OH36NFNTU6TTabxeL5ZlobXmzJkzVFZW0tHRgWHILQohxOKRQL+G4gzcMAzcbrf9UVNTg9/vZ8+ePeRyOXsGr5Qik8nQ09OD0+mkra1NQl0IsWgk0K+hWCPXWqOUQmuNZVm43W7a2trIZDLs3buXfD5POp0GIJlMYpomBw4cAJBQF0IsGgn0a3C73YRCIcLhsF0rD4VC9pL/jo4OkskkPT09dvB7vV601kQiEQ4dOgRIqAshFocE+jUopaiursbtdpNIJPD5fASDQfsGqWmabNu2DaUU3d3dAPZeL9lsFq01hw4dwul00tLSIu2MQohbSqaN16C1ZnJykrGxMeLxOGNjY0xOTl60gMg0TTZs2MCaNWvsEk0+n7fDe3p6miNHjsjCIyHELSeBfg3pdJpwOIzf78fn8+H3+wmHw3a9vMjr9dLR0UFNTQ1KKSzLwuFwkM/nSSaTDA4Osm/fPnvFqRBC3AoS6NdQ7HIpzraVUhiGQTabvehxSinq6urYuXMnzc3N+Hw+HA4HuVzO7n4ZHx/n5MmTJJPJpbgUIcRtQAL9Gi7scgHsLhen03nZY5VShEIhOjs7CYVCAPZMXWtNJpOhv7+fs2fPYlnWol6HEOL2IDdFr+FaXS5XopSiubmZbDbLoUOHmJ6etvd9SSQSaK05cuQIXq9XOl+EEAvuJQNdKfV54PXAuNa6c+5zlcA3gFagD/gtrfXMrRvm0gkGg3bw+v1+PB7PNbtVDMOgra0NgL179zI9PY1lWfYCpXg8Lu2MQohb4nrS5IvAay753IeBp7TW64Cn5v67pGitmZiYYHBwkHA4zPT0NLFY7Lr+32Kod3Z2EggE8Pv99otAMplkZmaGAwcO0NvbK+UXIcSCeclA11r/DJi+5NNvAL409/svAW9c4HEtuevtcLkawzBob2+ntbUVl8tlB7rD4cCyLHvhkYS6EGKh3GwNvU5rPQKgtR5RStUu4JiWhWt1uHg8nuv6GsV2xnQ6TTweJ5/P2+GdSqVk4ZEQYkHd8puiSqmHgYcBWlpabvXTzYtlWcQiEVLnz6Oqq8nn85ft43KlDperubCdUSnFyMgI8XjcbmeEFxcelZWVydmkQoh5udlAH1NKNczNzhuA8as9UGv9KPAoQFdX17JdVWNZFn0//CGr3vUuysfHsZxOgvX1WC0t5JqayDQ1EVi3Dvcdd0BbG9TUwHWEb7GdccuWLWit7RJL8TzSTCZjLzySE4+EEPNxs4H+feAh4JG5X7+3YCNaIomf/Yzmt70NXC6G3/c+zOlp9LlzVEQieJ98EmNq6uL/weeD1tbCx5o1l/9aUWEH/oUz9Ww2y/DwMFpr0un0ZQuPysrK8Pl8i3z14qUUf17ZbNY+vUpeeMVycz1ti48DrwCqlVJDwJ9RCPJvKqXeCQwAv3krB3nLPfMMvl/7NdJ+Pyf/8R/JrV5NMplkcnKSxsZGQqEQFQ4H1fE4qq8P+vrg3LkXf/3FLyASufhrBoMXBbxqbSXU2sp2j4eMy8VEOo1lWRcdPN3X10dZWRlbt26VE4+WkWLH06XrEYpbPQixXLxkoGut33yVP3r1Ao/lltJak0wmmZ6eJh6PY1kW2WyW4HPPsfp97yNVU8O/ve99OJTC6usjmUyilMLhcODz+ZhJJAiuW4dny5YrP0E4fHHIF3/t7YWnnoJ4HAWsAn4dSAcCRCsriVVXM1tVxUwoRKyqiv7z5yEeZ+vdd0uo3wI3M9NOpVKMj4/jdrtxOBw4HA7C4TBlZWXXfYNciMVwW6wU1VozNjZGT08PIyMjhMNhkskka48c4T89/jgTtbV886GHSFsW+uxZAPL5PMFgkOPHjxOLxaivr792h0soBNu3Fz4uHwBMTdlBb/X2kuruJnfsGJXDw7QcOYLjkv1hclVV6PZ21JXKOatXw1VWq4qLaa1JpVLE43H7Rb14CMn1zLS11oyOjjI5OWnvdR8MBnG5XDfU8STEYii5QL90BuZyuYjFYpw7d45wOIxpmhiGweaDB3ngm99kpKmJ//2Od5D2+fC53eRyOfL5vD17y2aznDt3DphHl45SUF1d+Ni5EwMIWhYTvb08c+AAkZkZnDMzlE1NEZiYIDg1Rc3sLKsyGbz79qG+8x24JPBZterKtfvWVmhuhhvoxllIy6nWrLVmfHycvr4+ZmdnSSaTZDIZ+4xYj8fzkjPtVCplb31cbDkdGRnB4XCQyWTsNQrBYPAlVxELcauVVKBfWuss9n3HYjFGR0eJRqNordn0zDO84pvfpL+tje8+9BCZCxb+FMssgN2HrrW292RZKMXVpNlslhdeeIFkKESyvJyJtWtJp9N4PB4qKyvZtm0bbatXY4yOXl7O6euDn/8cHn8cLlycZBjQ1HTlwF+zpvBicAvKOZZlcf78eSYnJ8lmsySTSTweDzU1NVRWVlJWVraoWx2k02nGxsaYnp4mm80yOztLOBwmHo+TSCSoqamxX7QvDfTizL63t5eJiQlSqRQDAwPkcjmi0SimadovVhUVFaxevZrVq1dTW1sroS6WTMkEutaaaDTK8PCw/dbYsixOnjyJy+UiHo8Tj8fZ9uSTvPz73+fsxo189y1vwXK54IJecygEk8fjoby8nGw2i8PhoKqq6paE+vr165menqanp8cOEcMwSKVS9hYBMLfvS3MzvPzll3+hbBaGhl4M+gtD/8knYXi4UPYpcjqhpeXqHTp1dYUXhRugteb8+fOcOnWKcDjMxMQEmUyGTCZDMBikoqKCtrY2tm3bZr9g3kqWZTE5OUlvb69dLgmHw2QyGftc2GQySXNz82VrC4oTg7GxMYaHh4lEIvbOmcUXh3w+j2ma+P1+EokE4XAYr9dLeXm5lGHEkimJQC/+AxwaGmJgYIB0Oo3f7yefzzM6OorX62U2FmPnD37A7iee4PjWrfzwwQfJGwaKwqzc6XTaJRoo1NDj8Thut5v6+nq8Xu8NLSq6XsVj7FKpFENDQxddk9aaWCzG0aNHcblcNDc3X3n253S+OPt+5Ssv//N0GgYGrnzT9gc/gLGxix/v8RTq9Fea3be2FkpHl4wjlUoxMjJCJBJhYmLCnqE7HA5isRj
|
||
|
|
"text/plain": [
|
||
|
|
"<Figure size 432x288 with 1 Axes>"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
"metadata": {
|
||
|
|
"needs_background": "light"
|
||
|
|
},
|
||
|
|
"output_type": "display_data"
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"name": "stdout",
|
||
|
|
"output_type": "stream",
|
||
|
|
"text": [
|
||
|
|
"[[Model]]\n",
|
||
|
|
" Model(piecewise_linear)\n",
|
||
|
|
"[[Fit Statistics]]\n",
|
||
|
|
" # fitting method = ampgo, with L-BFGS-B as local solver\n",
|
||
|
|
" # function evals = 79866\n",
|
||
|
|
" # data points = 491\n",
|
||
|
|
" # variables = 8\n",
|
||
|
|
" chi-square = 65.4379204\n",
|
||
|
|
" reduced chi-square = 0.13548224\n",
|
||
|
|
" Akaike info crit = -973.533027\n",
|
||
|
|
" Bayesian info crit = -939.961474\n",
|
||
|
|
"[[Variables]]\n",
|
||
|
|
" x0: 62.1499306 (init = 40.5)\n",
|
||
|
|
" x1: 181.207436 (init = 40.5)\n",
|
||
|
|
" x2: 196.301872 (init = 40.5)\n",
|
||
|
|
" b: 5.81709441 (init = 0)\n",
|
||
|
|
" k1: 0.06544372 (init = 0)\n",
|
||
|
|
" k2: -0.07421739 (init = 0)\n",
|
||
|
|
" k3: -0.32199700 (init = 0)\n",
|
||
|
|
" k4: 0.29325377 (init = 0)\n",
|
||
|
|
"\n"
|
||
|
|
]
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"source": [
|
||
|
|
"from lmfit import Model, Parameters\n",
|
||
|
|
"\n",
|
||
|
|
"def get_data():\n",
|
||
|
|
" site_id='NARRA0018'\n",
|
||
|
|
" profile_type='prestorm'\n",
|
||
|
|
" df_profile = df_profiles.query(\"site_id == '{}' and profile_type == '{}'\".format(site_id, profile_type))\n",
|
||
|
|
" x = np.array(df_profile.index.get_level_values('x').tolist())\n",
|
||
|
|
" z = np.array(df_profile.z.tolist()) \n",
|
||
|
|
"\n",
|
||
|
|
" nan_mask = ma.masked_invalid(z).mask\n",
|
||
|
|
" x = x[~nan_mask]\n",
|
||
|
|
" z = z[~nan_mask]\n",
|
||
|
|
" return x,z\n",
|
||
|
|
"\n",
|
||
|
|
"# def piecewise_linear(x, x0, x1, b, k1, k2, k3):\n",
|
||
|
|
"# condlist = [x < x0, (x >= x0) & (x < x1), x >= x1]\n",
|
||
|
|
"# funclist = [lambda x: k1*x + b, lambda x: k1*x + b + k2*(x-x0), lambda x: k1*x + b + k2*(x-x0) + k3*(x - x1)]\n",
|
||
|
|
"# return np.piecewise(x, condlist, funclist)\n",
|
||
|
|
"\n",
|
||
|
|
"# x,z = get_data()\n",
|
||
|
|
"\n",
|
||
|
|
"# fmodel = Model(piecewise_linear)\n",
|
||
|
|
"# params = Parameters()\n",
|
||
|
|
"# params.add('x0', value=0, vary=True, min=min(x), max=max(x))\n",
|
||
|
|
"# params.add('x1', value=0, vary=True, min=min(x), max=max(x))\n",
|
||
|
|
"# params.add('b', value=0, vary=True)\n",
|
||
|
|
"# params.add('k1', value=0, vary=True, min=-0.01, max=0.01)\n",
|
||
|
|
"# params.add('k2', value=0, vary=True, min=-0.1, max=-0.5)\n",
|
||
|
|
"# params.add('k3', value=0, vary=True, min=0.1, max=0.5)\n",
|
||
|
|
"\n",
|
||
|
|
"def piecewise_linear(x, x0, x1, x2, b, k1, k2, k3,k4):\n",
|
||
|
|
" condlist = [x < x0, (x >= x0) & (x < x1), (x >= x1) & (x < x2), x >= x2]\n",
|
||
|
|
" funclist = [lambda x: k1*x + b, lambda x: k1*x + b + k2*(x-x0), lambda x: k1*x + b + k2*(x-x0) + k3*(x - x1), lambda x: k1*x + b + k2*(x-x0) + k3*(x - x1) +k4*(x-x2)]\n",
|
||
|
|
" return np.piecewise(x, condlist, funclist)\n",
|
||
|
|
"\n",
|
||
|
|
"x,z = get_data()\n",
|
||
|
|
"\n",
|
||
|
|
"fmodel = Model(piecewise_linear)\n",
|
||
|
|
"params = Parameters()\n",
|
||
|
|
"params.add('x0', value=0, vary=True, min=min(x), max=max(x))\n",
|
||
|
|
"params.add('x1', value=0, vary=True, min=min(x), max=max(x))\n",
|
||
|
|
"params.add('x2', value=0, vary=True, min=min(x), max=max(x))\n",
|
||
|
|
"params.add('b', value=0, vary=True)\n",
|
||
|
|
"params.add('k1', value=0, vary=True, min=-0.5, max=0.5)\n",
|
||
|
|
"params.add('k2', value=0, vary=True, min=-0.5, max=0.5)\n",
|
||
|
|
"params.add('k3', value=0, vary=True, min=-0.5, max=0.5)\n",
|
||
|
|
"params.add('k4', value=0, vary=True, min=-0.5, max=0.5)\n",
|
||
|
|
"\n",
|
||
|
|
"\n",
|
||
|
|
"result = fmodel.fit(z, params, x=x,method='ampgo')\n",
|
||
|
|
"\n",
|
||
|
|
"\n",
|
||
|
|
"pyplot.figure()\n",
|
||
|
|
"pyplot.plot(x,z, color='0.5',marker='.', alpha=.2, ms=10,linestyle=\"None\")\n",
|
||
|
|
"pyplot.plot(x,result.best_fit, color='r')\n",
|
||
|
|
"pyplot.show()\n",
|
||
|
|
"print(result.fit_report())"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"metadata": {},
|
||
|
|
"source": [
|
||
|
|
"## Try spline"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 65,
|
||
|
|
"metadata": {
|
||
|
|
"ExecuteTime": {
|
||
|
|
"end_time": "2018-12-05T01:49:30.968871Z",
|
||
|
|
"start_time": "2018-12-05T01:49:30.648994Z"
|
||
|
|
},
|
||
|
|
"code_folding": [
|
||
|
|
2
|
||
|
|
]
|
||
|
|
},
|
||
|
|
"outputs": [
|
||
|
|
{
|
||
|
|
"data": {
|
||
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD8CAYAAABn919SAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzs3Xl01Nd9///nnX3RMtr3HRAgAcJIBgyxHcdL4t2O46xO7KZJ47iJ03N6epqTnvbb7zenJ9/fN2lip2lax67TxIntxguO45J4i+MFGyNAgFgESGjfNZrRaDT73N8faKYYGwxIaBm9H+foIA0f6XM/GumlO/fzvvcqrTVCCCEWP8N8N0AIIcTskEAXQogUIYEuhBApQgJdCCFShAS6EEKkCAl0IYRIERLoQgiRIiTQhRAiRUigCyFEijDN5clyc3N1ZWXlXJ5SCCEWvd27d49qrfM+7Lg5DfTKykqam5vn8pRCCLHoKaW6zuU4GXIRQogUIYEuhBApQgJdCCFShAS6EEKkCAl0IYRIEXNa5bLYaa0JhUJEIhHMZjNWqxWl1Hw3SwghAAn0cxaPx+nt7WV4eBgAq9VKfn4++fn5EupCiAVBAv0cxONx2tvbaWlpIR6PYzQaMZvNjIyMYLFYcLlcEupCiHkngf4htNb09vayb98+xsbGUEoRiURQSmE2m4nFYqxZs4aCggIJdSHEvJJA/xChUIihoSG01lgsFiKRCJFIBGM4TNWxY8S6utgfj9N06aXSUxdCzCsJ9A8RiUQwGo2YTCaMRiOBQAA1OckdDz9McXc3AIffeIPmv/s71qxdKz11IcS8kUD/ECbTyW9RoqrFPzHBTb/+NYW9vbxy113k9PXR8OqrjDz7LAeNRtLT03E6nfPcaiHEUiR16GehtWZiYiI5Zg7QtGMHNUeO8PqnPkXrhg28efPNDJaXc8mTT3Li0CFaWlqIxWLz3HIhxFIkgX4WoVAIr9dLSUkJq1evZn1hIZf89rdMbN3K6Cc/ic1mI64Uuz71KdLGx6l/9VWOHz/O0aNHicfj8918IcQSI4F+FpFIBIPBgMFgwGq1UvrIIxgiEdSDD9Kwfj1paWkYjUa6q6o4sWYN63//ewxjY7S2ttLR0SGhLoSYUxLoZ2E2m4nH42itMYyN4fjNb/DefjvmVauorq5m/XSo2+123rn1VszhMBsfewzfxAR79uyRUBdCzCkJ9LOwWq24XC78fj/mhx9GhULEvvENrFYrBoOB2tpa6uvrMRqNuAsK2Hn99azYu5dlb76J1+ulpaVFQl0IMWekyuUslFLk5uZiNRpJe+IJItddR9bmzckbpEajkYaGBpRS7N+/n6O33Ub5sWNc9thjDLhcuFetoqWlBbPZTHl5uZQzCiEuKumhn4XWmtHRUSb/678wDg8zfMstjI6OorVOHmM0Glm5ciVVVVWYrFb+8OUvM+lyceMjj+B0u3G73Rw4cACPx/OezxNCiNkmgX4WoVAIj8dD9rZtxPLz4frr8Xg8hEKh9xxnt9upra0lLy+PcHo6v/3zP8cUjfKJBx8kNjZGT08Pzc3NyRmnQghxMUign0UkEsE8OortlVeYuuMOlNmMwWAgEom85zilFAUFBTQ1NVFWVkaoupoX772XzOFhbvqP/8AYjTI8PExbWxuBQGCerkYIkeok0M/CbDaT9swzqFiMqU9/Gq018Xgcs9n8vmOVUrhcLurr63G5XAzU1vLi5z5H6fHjXPXznxMOBunq6qK9vV1ukgohLgoJ9LOwWixkPfccgUsuwVdcjN/vx+VyYbVaP/B4pRRlZWWsX7+e7OxserZuZcdtt7Fszx7W//rXTExMcODAAal8EUJcFB9a5aKU+g/gRmBYa10//Vg28CRQCXQCd2qtxy9eM+fJzp2Yjh8n8MMf4nK5cDqd2Gy2s1arGAwGqqurAdi1axd7r7oKp9tNw6uv4s/K4tB119HS0gJAdXU1BoP8TRVCzI5zSZOfAx8/7bG/BV7RWi8HXpn+OKVorQn8+78Tt9kYvuIK3G43Pp/vnD43Eer19fWkpaez9667aF+/ni1PP03pjh2Mj4/LxCMhxKz70EDXWr8OuE97+BbgP6ff/0/g1llu17wLeb1Yn3mG4PXXYy8owOl0fmCFy5kYDAZqamqorKzEbLPx6j33MFBTw7WPPUbBkSMy8UgIMesu9PV+gdZ6AGD63/zZa9LCoLdtw+jzMXXnncDJ8fEPqnA5m1PLGbXVyvavfQ1vbi7X/vSn2NrbcbvdtLS00NPTI+WMQogZu+gDuEqpryqlmpVSzSMjIxf7dDMSj8fxer0MDAxg+OUviRQVEbrsMoCzVricyanljFVVVZjy83nhvvuIms3c+u//jmN8XCYeCSFmzYUG+pBSqghg+t/hMx2otX5Ia92otW7My8u7wNNdfPF4nKNHj7Jnzx6633oLy2uvMfzxjzM5NcXU1NSHVricSaKccc2aNRQWFuLPzeW5r34VazB4cuKR2y0Tj4QQs+JCA/23wJem3/8S8NzsNGf++Hw+BgYGyMrKouqFFwBov/LK5AJdZWVl5OXlXdB6LKdPPJpasYKXvvY1XIOD3CgTjxYFrTXBYBCfz0cwGJQ/vGJB+tBAV0o9DrwN1CqlepVSXwa+B1yjlDoGXDP98aLm9/uJxWJExsfJfvppRq+4grH0dIaHh8+rwuVMTp941LdqFS997nOUHTvGR3/xCyKhEJ2dnbS1tcmORwuM1pqRkRG6u7sZGhqiu7ubkZERCXWx4HxoHbrW+rNn+K+PzXJbLiqtNYFAALfbjd/vJx6PJ7eWMxqNDA8P09vbS+HLL2Py+WjeuhW/34/JZMLhcODxeMjIyMBms11wGxITjyKRyMmboZdfzk6/n43PPovv8cd5+9Zb2b9/P/F4nIaGBoxG4yx+BwSc/DkIhUInl3WY3if2w151BYNBhoeHsVqtmEwmTCbTrPw8CDHblsTyuVprhoaGaG1tZWBgAI/HQyAQQGuNUirZ0yoaGKD2ySc5Wl/PPoeDdK+Xw4cP4/P5KCwsJBKJzPgX+NSJR3v27KH1+uuxj4xwycsv48/K4vgnPkFrayt2u52VK1fKxKMZSgyV+P3+5B/1QCCA0WgkHo/jcrnOOpSmtWZwcJDR0VHsdjtaa9LT07FYLLPy8yDEbEq5QD+9B2axWPD5fJw4cQKPx4PRaMSoNXV791LQ3U3YZmN/YyNxg4Fbf/ELAmlpvPzJTyZ7b5FIhBMnTgBQXl4+K208PdTf/cIXSPP52PrUU0xmZNDe0MDevXuxWCyLbjbphfSAL2ZbhoeH6ezsZHJykkAgQDgcpqSkhPz8fGw224f2tIPBYLICKTFfYGBgAJPJRDgcxuFw4HQ6SU9P/9BZxEJcbCkV6ImxTo/Hg8FgIBaLEY/H8fl8jLa3E+vvJ3NkhGsef5y8vj7CFgumaJTNf/wjcaWIWiw8fe+9hF0ujJwM3kQPPhqNzmpbE6EeiUR49913eeXP/oybfvQjrvnP/2TS6WS8rm7RLREQj8fp6+tjdHSUSCRCIBDAZrORl5dHdnY2GRkZc3odoVCIoaEh3G43kUiEyclJPB4Pfr+fqakp8vLykn+0Tw/0RM++o6ODkZERgsEg3V1dFOzfj7GvjxOrVtFWVIRSiqysLCoqKqioqCA/P19CXcyblAl0rTUTExP09/cnXxrH43GOtbay5te/ZsPzz2Ocvtnoc7n47Re+wOG6OtInJ1m7cyeWSITDW7cynptLPB7HZrORmZlJJBLBZDKRk5NzUUJ9xYoVuN1uWltb+e977+W2//f/uOnhh/nN/fczrhR79uwBFn6oa63p6+vj6NGjeDweRkZGCIfDhMNh0tPTycrKorq6moaGBkymi/9jF4/HGR0dpaOjIzlc4vF4CIfDxONxrFYrgUCAsrKy980tSHQMhoaG6O/vx+v14uzt5cZ/+Rdypl+tRcxm3vrkJ+m89lqmpqbweDzY7XYyMzNlGEbMm5QI9MQvYG9vL93d3YRCIZxOJwaPh8Z//EeKjx/n0ObN9JWUENCa4+vWEbV
|
||
|
|
"text/plain": [
|
||
|
|
"<Figure size 432x288 with 1 Axes>"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
"metadata": {
|
||
|
|
"needs_background": "light"
|
||
|
|
},
|
||
|
|
"output_type": "display_data"
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"source": [
|
||
|
|
"from scipy.signal import savgol_filter\n",
|
||
|
|
"\n",
|
||
|
|
"def get_data():\n",
|
||
|
|
" site_id='NARRA0018'\n",
|
||
|
|
" profile_type='prestorm'\n",
|
||
|
|
" df_profile = df_profiles.query(\"site_id == '{}' and profile_type == '{}'\".format(site_id, profile_type))\n",
|
||
|
|
" x = np.array(df_profile.index.get_level_values('x').tolist())\n",
|
||
|
|
" z = np.array(df_profile.z.tolist()) \n",
|
||
|
|
"\n",
|
||
|
|
" nan_mask = ma.masked_invalid(z).mask\n",
|
||
|
|
" x = x[~nan_mask]\n",
|
||
|
|
" z = z[~nan_mask]\n",
|
||
|
|
" return x,z\n",
|
||
|
|
"\n",
|
||
|
|
"x,z = get_data()\n",
|
||
|
|
"\n",
|
||
|
|
"z_filtered = savgol_filter(z, 31, 3)\n",
|
||
|
|
"\n",
|
||
|
|
"\n",
|
||
|
|
"pyplot.figure()\n",
|
||
|
|
"pyplot.plot(x,z, color='0.5',marker='.', alpha=.2, ms=10,linestyle=\"None\")\n",
|
||
|
|
"pyplot.plot(x,z_filtered, color='r')\n",
|
||
|
|
"pyplot.show()\n"
|
||
|
|
]
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"metadata": {
|
||
|
|
"kernelspec": {
|
||
|
|
"display_name": "Python 3",
|
||
|
|
"language": "python",
|
||
|
|
"name": "python3"
|
||
|
|
},
|
||
|
|
"language_info": {
|
||
|
|
"codemirror_mode": {
|
||
|
|
"name": "ipython",
|
||
|
|
"version": 3
|
||
|
|
},
|
||
|
|
"file_extension": ".py",
|
||
|
|
"mimetype": "text/x-python",
|
||
|
|
"name": "python",
|
||
|
|
"nbconvert_exporter": "python",
|
||
|
|
"pygments_lexer": "ipython3",
|
||
|
|
"version": "3.6.7"
|
||
|
|
},
|
||
|
|
"toc": {
|
||
|
|
"base_numbering": 1,
|
||
|
|
"nav_menu": {},
|
||
|
|
"number_sections": true,
|
||
|
|
"sideBar": true,
|
||
|
|
"skip_h1_title": false,
|
||
|
|
"title_cell": "Table of Contents",
|
||
|
|
"title_sidebar": "Contents",
|
||
|
|
"toc_cell": false,
|
||
|
|
"toc_position": {},
|
||
|
|
"toc_section_display": true,
|
||
|
|
"toc_window_display": false
|
||
|
|
},
|
||
|
|
"varInspector": {
|
||
|
|
"cols": {
|
||
|
|
"lenName": 16,
|
||
|
|
"lenType": 16,
|
||
|
|
"lenVar": 40
|
||
|
|
},
|
||
|
|
"kernels_config": {
|
||
|
|
"python": {
|
||
|
|
"delete_cmd_postfix": "",
|
||
|
|
"delete_cmd_prefix": "del ",
|
||
|
|
"library": "var_list.py",
|
||
|
|
"varRefreshCmd": "print(var_dic_list())"
|
||
|
|
},
|
||
|
|
"r": {
|
||
|
|
"delete_cmd_postfix": ") ",
|
||
|
|
"delete_cmd_prefix": "rm(",
|
||
|
|
"library": "var_list.r",
|
||
|
|
"varRefreshCmd": "cat(var_dic_list()) "
|
||
|
|
}
|
||
|
|
},
|
||
|
|
"types_to_exclude": [
|
||
|
|
"module",
|
||
|
|
"function",
|
||
|
|
"builtin_function_or_method",
|
||
|
|
"instance",
|
||
|
|
"_Feature"
|
||
|
|
],
|
||
|
|
"window_display": false
|
||
|
|
}
|
||
|
|
},
|
||
|
|
"nbformat": 4,
|
||
|
|
"nbformat_minor": 2
|
||
|
|
}
|