{ "cells": [ { "cell_type": "code", "execution_count": 1, "id": "22780974", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import os\n", "import matplotlib.pylab as plt\n", "import matplotlib\n", "import cartopy.crs as ccrs\n", "import cmocean.cm as cm\n", "import cartopy.feature as cfeature\n", "\n", "from matplotlib.animation import FuncAnimation\n", "import matplotlib.animation as animation\n", "from IPython.display import Video" ] }, { "cell_type": "code", "execution_count": 2, "id": "07157b22", "metadata": {}, "outputs": [], "source": [ "# font size for plotting\n", "fs = 30\n", "font = {'size' : 18}\n", "matplotlib.rc('font', **font)\n", "\n", "# extend of the grid:\n", "minlat = 20\n", "maxlat = 50\n", "minlon = 100\n", "maxlon = 200\n", "# define the grid\n", "X,Y = np.meshgrid(np.arange(minlon,maxlon,3), np.arange(minlat, maxlat, 2))\n", "# define time and Earth rotation speed\n", "rotateperday = 0.5 # degrees per day\n", "days = np.arange(0,180,0.15)\n", "midlons = 100 + rotateperday*days\n", "\n", "filedir = os.path.abspath(\".\")\n", "vdir = os.path.abspath(os.path.join(filedir, \"..\", \"images\"))\n", "vfile = 'RotatingEarth.mp4'\n", "vfile_path = os.path.join(vdir, vfile)" ] }, { "cell_type": "code", "execution_count": 3, "id": "c2f66989", "metadata": {}, "outputs": [], "source": [ "# Functions that define the time varying flow field\n", "def xcon(x, lat, r0):\n", " x = (x - minlon) / (maxlon-minlon) * r0*np.pi\n", " return x\n", "\n", "def ycon(y):\n", " y = (y-minlat) / (maxlat-minlat) * 8000 - 4000\n", " return y\n", "\n", "def BickleyJetflow(x, y, t, midlat = 35):\n", " \"\"\"\n", " Return the Bickley Jet flow field\n", "\n", " \"\"\"\n", " ua = np.zeros(x.shape)\n", " va = np.zeros(x.shape)\n", " \n", " Ubj = 0.06266\n", " L = 1700.\n", " r0 = 6371.\n", " k1 = 2 * 1/ r0\n", " k2 = 2 * 2 / r0\n", " k3 = 2 * 3/ r0\n", " eps1 = 0.075\n", " eps2 = 0.4\n", " eps3 = 0.3\n", " c3 = 0.461 * Ubj\n", " c2 = 0.205 * Ubj\n", " c1 = 0.1446 * Ubj\n", " \n", " x = xcon(x, midlat, r0=r0)\n", " y = ycon(y)\n", " \n", " def v(x1, x2, t): #2D Bickley Jet velocity given single location\n", " f1 = eps1 * np.exp(-1j *k1 * c1 * t)\n", " f2 = eps2 * np.exp(-1j *k2 * c2 * t)\n", " f3 = eps3 * np.exp(-1j *k3 * c3 * t)\n", " F1 = f1 * np.exp(1j * k1 * x1)\n", " F2 = f2 * np.exp(1j * k2 * x1)\n", " F3 = f3 * np.exp(1j * k3 * x1)\n", " G = np.real(np.sum([F1,F2,F3]))\n", " G_x = np.real(np.sum([1j * k1 *F1, 1j * k2 * F2, 1j * k3 * F3]))\n", " u = Ubj / (np.cosh(x2/L)**2) + 2 * Ubj * np.sinh(x2/L) / (np.cosh(x2/L)**3) * G\n", " vd = Ubj * L * (1./np.cosh(x2/L))**2 * G_x\n", " return [u,vd]\n", " \n", " for i in range(x.shape[0]):\n", " for j in range(y.shape[1]):\n", " ua[i,j], va[i,j] = v(x[i,j], y[i,j], t*36000)\n", " res = np.sqrt(ua**2+va**2)\n", "\n", " return res" ] }, { "cell_type": "code", "execution_count": 4, "id": "e6c01774", "metadata": {}, "outputs": [], "source": [ "# The initial plot\n", "\n", "def plot(t): \n", " data = BickleyJetflow(X, Y, days[t])\n", " \n", " fig = plt.figure(figsize=(6,5))\n", " projection = ccrs.NearsidePerspective(midlons[t], 20)\n", "\n", " ax = plt.subplot(111, projection=projection)\n", " ax.add_feature(cfeature.OCEAN, zorder=0)\n", " ax.add_feature(cfeature.LAND, zorder=0, edgecolor='black')\n", " ax.set_global()\n", " tit = ax.set_title('day %.1f'%(days[t]), fontsize=fs)\n", "\n", " im = ax.pcolor(X, Y, data, transform=ccrs.PlateCarree(), cmap=cm.speed,\n", " vmin=0, vmax=0.1)\n", "\n", " #% add colorbar \n", " fig.subplots_adjust(right=0.8)\n", " cbar_ax = fig.add_axes([0.8, 0.15, 0.1, 0.7])\n", " cbar_ax.set_visible(False)\n", " \n", " cbar = fig.colorbar(im, ax=cbar_ax, orientation = 'vertical', fraction = 0.8,\n", " aspect=18, extend='max')\n", " cbar.ax.xaxis.set_label_position('bottom')\n", " cbar.ax.set_xlabel('m/s', fontsize=fs)\n", " return fig, tit, cbar_ax" ] }, { "cell_type": "code", "execution_count": 5, "id": "c7e5fe76", "metadata": {}, "outputs": [], "source": [ "# The animation\n", "\n", "def animate(t):\n", " if(t%25==0):\n", " print(t/len(days),'%')\n", " data = BickleyJetflow(X, Y, days[t])\n", " projection = ccrs.NearsidePerspective(midlons[t], 20)\n", "\n", " ax = plt.subplot(111, projection=projection)\n", " ax.add_feature(cfeature.OCEAN, zorder=0)\n", " ax.add_feature(cfeature.LAND, zorder=0, edgecolor='black')\n", " ax.set_global()\n", " tit = ax.set_title('day %.1f'%(days[t]), fontsize=fs)\n", "\n", " ax.set_transform(projection)\n", " tit.set_text('day %.1f'%(days[t]))\n", " im = ax.pcolor(X, Y, data, transform=ccrs.PlateCarree(), cmap=cm.speed,\n", " vmin=0, vmax=0.1)\n", "\n", "if(__name__=='__main__'):\n", " if not os.path.exists(vfile_path):\n", " fig, tit, cbar_ax = plot(0)\n", " anim = FuncAnimation(fig, animate, frames = np.arange(len(days)), interval=10)\n", " Writer = animation.writers['ffmpeg']\n", " writer = Writer(fps=20, metadata=dict(artist='Me'), bitrate=1250)\n", " anim.save(vfile_path, writer=writer, dpi=92)\n", " plt.close()" ] }, { "cell_type": "code", "execution_count": 6, "id": "38122a11", "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Video(vfile_path, embed=True)" ] }, { "cell_type": "code", "execution_count": null, "id": "4ae4fc75", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python [conda env:Cart2SNS]", "language": "python", "name": "conda-env-Cart2SNS-py" }, "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.7.10" } }, "nbformat": 4, "nbformat_minor": 5 }