{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Analysing WRF output example\n", "We will use the wrf output for one day, namely 07-01-2015. We will read the data, check which variables are in the file, try to make some time series, and plot the data." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Reading the data and make timeseries plot" ] }, { "cell_type": "code", "execution_count": 41, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Times\n", "U\n", "V\n", "W\n", "PH\n", "PHB\n", "T\n", "P\n", "PB\n", "P_HYD\n", "Q2\n", "T2\n", "TH2\n", "PSFC\n", "U10\n", "V10\n", "XTIME\n", "QVAPOR\n", "QCLOUD\n", "QRAIN\n", "TSLB\n", "SMOIS\n", "GRDFLX\n", "SNOW\n", "SNOWH\n", "SSTSK\n", "QKE\n", "TKE_PBL\n", "EL_PBL\n", "MAPFAC_UX\n", "MAPFAC_UY\n", "HGT\n", "TSK\n", "RAINNC\n", "SNOWNC\n", "GRAUPELNC\n", "HAILNC\n", "CLDFRA\n", "SWDOWN\n", "GLW\n", "SWNORM\n", "OLR\n", "ALBEDO\n", "EMISS\n", "NOAHRES\n", "UST\n", "PBLH\n", "HFX\n", "QFX\n", "LH\n", "SNOWC\n", "SR\n", "SST\n" ] } ], "source": [ "# Import packages\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "import pandas as pd\n", "from netCDF4 import Dataset,num2date\n", "\n", "#Don't use this command in the terminal!!!\n", "%matplotlib notebook \n", "\n", "\n", "# Make the file object\n", "wrfout = Dataset('wrfout_d03_2015-07-01_00:00:00')\n", "\n", "# Show which variables are in the file\n", "for x in wrfout.variables:\n", " print x" ] }, { "cell_type": "code", "execution_count": 42, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "float32 XTIME(Time)\n", " FieldType: 104\n", " MemoryOrder: 0 \n", " description: minutes since 2015-06-25 00:00:00\n", " units: minutes since 2015-06-25 00:00:00\n", " stagger: \n", "unlimited dimensions: Time\n", "current shape = (24,)\n", "filling off" ] }, "execution_count": 42, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Inspect time array in more detail\n", "wrfout.variables['XTIME']" ] }, { "cell_type": "code", "execution_count": 43, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[ 8640. 8700. 8760. 8820. 8880. 8940. 9000. 9060. 9120.\n", " 9180. 9240. 9300. 9360. 9420. 9480. 9540. 9600. 9660.\n", " 9720. 9780. 9840. 9900. 9960. 10020.]\n" ] } ], "source": [ "# Extract time data from wrfout object\n", "wrftime = wrfout.variables['XTIME'][:]\n", "print wrftime" ] }, { "cell_type": "code", "execution_count": 44, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[datetime.datetime(2015, 7, 1, 0, 0) datetime.datetime(2015, 7, 1, 1, 0)\n", " datetime.datetime(2015, 7, 1, 2, 0) datetime.datetime(2015, 7, 1, 3, 0)\n", " datetime.datetime(2015, 7, 1, 4, 0) datetime.datetime(2015, 7, 1, 5, 0)\n", " datetime.datetime(2015, 7, 1, 6, 0) datetime.datetime(2015, 7, 1, 7, 0)\n", " datetime.datetime(2015, 7, 1, 8, 0) datetime.datetime(2015, 7, 1, 9, 0)\n", " datetime.datetime(2015, 7, 1, 10, 0) datetime.datetime(2015, 7, 1, 11, 0)\n", " datetime.datetime(2015, 7, 1, 12, 0) datetime.datetime(2015, 7, 1, 13, 0)\n", " datetime.datetime(2015, 7, 1, 14, 0) datetime.datetime(2015, 7, 1, 15, 0)\n", " datetime.datetime(2015, 7, 1, 16, 0) datetime.datetime(2015, 7, 1, 17, 0)\n", " datetime.datetime(2015, 7, 1, 18, 0) datetime.datetime(2015, 7, 1, 19, 0)\n", " datetime.datetime(2015, 7, 1, 20, 0) datetime.datetime(2015, 7, 1, 21, 0)\n", " datetime.datetime(2015, 7, 1, 22, 0) datetime.datetime(2015, 7, 1, 23, 0)]\n" ] } ], "source": [ "# Convert time to datetime format:\n", "timeunits = wrfout.variables['XTIME'].units\n", "time = num2date(wrftime,units=timeunits,calendar='standard')\n", "print time" ] }, { "cell_type": "code", "execution_count": 45, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "float32 T(Time, bottom_top, south_north, west_east)\n", " FieldType: 104\n", " MemoryOrder: XYZ\n", " description: perturbation potential temperature (theta-t0)\n", " units: K\n", " stagger: \n", " coordinates: XLONG XLAT XTIME\n", "unlimited dimensions: Time\n", "current shape = (24, 39, 201, 201)\n", "filling off\n", "\n" ] } ], "source": [ "# Import other variables, e.g. Temperature\n", "# let's inspect the variable first\n", "print wrfout.variables['T']" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can see that temperature has 4 dimensions. Also, the order in which the data is stored along the dimensions is given. We could just import the full matrix and check that it has indeed four dimensions as follows:" ] }, { "cell_type": "code", "execution_count": 46, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "(24, 39, 201, 201)" ] }, "execution_count": 46, "metadata": {}, "output_type": "execute_result" } ], "source": [ "wrfT = wrfout.variables['T'][:]\n", "wrfT.shape" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, since we are interested in plotting a time series, we need only one point in space, so we need to limit x,y and z, while keeping all the time units. Let's assume, for now, that our point of interest is x=100, y=100, z=0. That is, a point in more or less the middle of the domain, near the surface. Note that we strictly follow the order of the dimensions as it is given in the variable description" ] }, { "cell_type": "code", "execution_count": 47, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(24,)\n", "[-13.25781155 -13.25263882 -13.24032211 -13.39459515 -13.47680187\n", " -13.59323025 -13.6786747 -13.72670937 -13.72151375 -13.48870277\n", " -13.28540421 -13.01741314 -12.58317089 -12.16792202 -11.94770718\n", " -11.79158592 -11.6240406 -11.64165592 -11.44226551 -11.11795044\n", " -10.88765812 -10.86215973 -10.82411957 -10.93723869]\n" ] } ], "source": [ "wrfT = wrfout.variables['T'][:,0,100,100]\n", "print wrfT.shape\n", "print wrfT" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It is important to try to minimize the amount of data that is actively loaded into memory. Therefore, it is way better to import only one temperature point than the complete 4-dimensional array. If you put to much data in memory, python may become very slow and even crash. Recall that the units of T were K, and now we find values of -13. This is clearly not possible. That is, because T is actually a perturbation temperature, and from the [wrf user guide](http://www2.mmm.ucar.edu/wrf/users/docs/user_guide_V3/users_guide_chap5.htm#special_fields), we need to add 300K to it to find actual temperature values. " ] }, { "cell_type": "code", "execution_count": 58, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[ 286.7421875 286.74737549 286.75967407 286.60540771 286.52319336\n", " 286.4067688 286.32131958 286.27328491 286.2784729 286.5112915\n", " 286.71459961 286.98257446 287.4168396 287.83209229 288.05230713\n", " 288.20840454 288.37594604 288.3583374 288.55773926 288.88204956\n", " 289.11233521 289.1378479 289.1758728 289.06277466]\n" ] } ], "source": [ "temp = wrfT+300.\n", "print temp" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now that we've successfully imported the temperature and time, let's make a plot of the time series" ] }, { "cell_type": "code", "execution_count": 59, "metadata": { "collapsed": false }, "outputs": [ { "data": { "application/javascript": [ "/* Put everything inside the global mpl namespace */\n", "window.mpl = {};\n", "\n", "mpl.get_websocket_type = function() {\n", " if (typeof(WebSocket) !== 'undefined') {\n", " return WebSocket;\n", " } else if (typeof(MozWebSocket) !== 'undefined') {\n", " return MozWebSocket;\n", " } else {\n", " alert('Your browser does not have WebSocket support.' +\n", " 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n", " 'Firefox 4 and 5 are also supported but you ' +\n", " 'have to enable WebSockets in about:config.');\n", " };\n", "}\n", "\n", "mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n", " this.id = figure_id;\n", "\n", " this.ws = websocket;\n", "\n", " this.supports_binary = (this.ws.binaryType != undefined);\n", "\n", " if (!this.supports_binary) {\n", " var warnings = document.getElementById(\"mpl-warnings\");\n", " if (warnings) {\n", " warnings.style.display = 'block';\n", " warnings.textContent = (\n", " \"This browser does not support binary websocket messages. \" +\n", " \"Performance may be slow.\");\n", " }\n", " }\n", "\n", " this.imageObj = new Image();\n", "\n", " this.context = undefined;\n", " this.message = undefined;\n", " this.canvas = undefined;\n", " this.rubberband_canvas = undefined;\n", " this.rubberband_context = undefined;\n", " this.format_dropdown = undefined;\n", "\n", " this.image_mode = 'full';\n", "\n", " this.root = $('
');\n", " this._root_extra_style(this.root)\n", " this.root.attr('style', 'display: inline-block');\n", "\n", " $(parent_element).append(this.root);\n", "\n", " this._init_header(this);\n", " this._init_canvas(this);\n", " this._init_toolbar(this);\n", "\n", " var fig = this;\n", "\n", " this.waiting = false;\n", "\n", " this.ws.onopen = function () {\n", " fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n", " fig.send_message(\"send_image_mode\", {});\n", " fig.send_message(\"refresh\", {});\n", " }\n", "\n", " this.imageObj.onload = function() {\n", " if (fig.image_mode == 'full') {\n", " // Full images could contain transparency (where diff images\n", " // almost always do), so we need to clear the canvas so that\n", " // there is no ghosting.\n", " fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n", " }\n", " fig.context.drawImage(fig.imageObj, 0, 0);\n", " };\n", "\n", " this.imageObj.onunload = function() {\n", " this.ws.close();\n", " }\n", "\n", " this.ws.onmessage = this._make_on_message_function(this);\n", "\n", " this.ondownload = ondownload;\n", "}\n", "\n", "mpl.figure.prototype._init_header = function() {\n", " var titlebar = $(\n", " '
');\n", " var titletext = $(\n", " '
');\n", " titlebar.append(titletext)\n", " this.root.append(titlebar);\n", " this.header = titletext[0];\n", "}\n", "\n", "\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "\n", "mpl.figure.prototype._root_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "mpl.figure.prototype._init_canvas = function() {\n", " var fig = this;\n", "\n", " var canvas_div = $('
');\n", "\n", " canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n", "\n", " function canvas_keyboard_event(event) {\n", " return fig.key_event(event, event['data']);\n", " }\n", "\n", " canvas_div.keydown('key_press', canvas_keyboard_event);\n", " canvas_div.keyup('key_release', canvas_keyboard_event);\n", " this.canvas_div = canvas_div\n", " this._canvas_extra_style(canvas_div)\n", " this.root.append(canvas_div);\n", "\n", " var canvas = $('');\n", " canvas.addClass('mpl-canvas');\n", " canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n", "\n", " this.canvas = canvas[0];\n", " this.context = canvas[0].getContext(\"2d\");\n", "\n", " var rubberband = $('');\n", " rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n", "\n", " var pass_mouse_events = true;\n", "\n", " canvas_div.resizable({\n", " start: function(event, ui) {\n", " pass_mouse_events = false;\n", " },\n", " resize: function(event, ui) {\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " stop: function(event, ui) {\n", " pass_mouse_events = true;\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " });\n", "\n", " function mouse_event_fn(event) {\n", " if (pass_mouse_events)\n", " return fig.mouse_event(event, event['data']);\n", " }\n", "\n", " rubberband.mousedown('button_press', mouse_event_fn);\n", " rubberband.mouseup('button_release', mouse_event_fn);\n", " // Throttle sequential mouse events to 1 every 20ms.\n", " rubberband.mousemove('motion_notify', mouse_event_fn);\n", "\n", " rubberband.mouseenter('figure_enter', mouse_event_fn);\n", " rubberband.mouseleave('figure_leave', mouse_event_fn);\n", "\n", " canvas_div.on(\"wheel\", function (event) {\n", " event = event.originalEvent;\n", " event['data'] = 'scroll'\n", " if (event.deltaY < 0) {\n", " event.step = 1;\n", " } else {\n", " event.step = -1;\n", " }\n", " mouse_event_fn(event);\n", " });\n", "\n", " canvas_div.append(canvas);\n", " canvas_div.append(rubberband);\n", "\n", " this.rubberband = rubberband;\n", " this.rubberband_canvas = rubberband[0];\n", " this.rubberband_context = rubberband[0].getContext(\"2d\");\n", " this.rubberband_context.strokeStyle = \"#000000\";\n", "\n", " this._resize_canvas = function(width, height) {\n", " // Keep the size of the canvas, canvas container, and rubber band\n", " // canvas in synch.\n", " canvas_div.css('width', width)\n", " canvas_div.css('height', height)\n", "\n", " canvas.attr('width', width);\n", " canvas.attr('height', height);\n", "\n", " rubberband.attr('width', width);\n", " rubberband.attr('height', height);\n", " }\n", "\n", " // Set the figure to an initial 600x600px, this will subsequently be updated\n", " // upon first draw.\n", " this._resize_canvas(600, 600);\n", "\n", " // Disable right mouse context menu.\n", " $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n", " return false;\n", " });\n", "\n", " function set_focus () {\n", " canvas.focus();\n", " canvas_div.focus();\n", " }\n", "\n", " window.setTimeout(set_focus, 100);\n", "}\n", "\n", "mpl.figure.prototype._init_toolbar = function() {\n", " var fig = this;\n", "\n", " var nav_element = $('
')\n", " nav_element.attr('style', 'width: 100%');\n", " this.root.append(nav_element);\n", "\n", " // Define a callback function for later on.\n", " function toolbar_event(event) {\n", " return fig.toolbar_button_onclick(event['data']);\n", " }\n", " function toolbar_mouse_event(event) {\n", " return fig.toolbar_button_onmouseover(event['data']);\n", " }\n", "\n", " for(var toolbar_ind in mpl.toolbar_items) {\n", " var name = mpl.toolbar_items[toolbar_ind][0];\n", " var tooltip = mpl.toolbar_items[toolbar_ind][1];\n", " var image = mpl.toolbar_items[toolbar_ind][2];\n", " var method_name = mpl.toolbar_items[toolbar_ind][3];\n", "\n", " if (!name) {\n", " // put a spacer in here.\n", " continue;\n", " }\n", " var button = $('');\n", " button.click(method_name, toolbar_event);\n", " button.mouseover(tooltip, toolbar_mouse_event);\n", " nav_element.append(button);\n", " }\n", "\n", " // Add the status bar.\n", " var status_bar = $('');\n", " nav_element.append(status_bar);\n", " this.message = status_bar[0];\n", "\n", " // Add the close button to the window.\n", " var buttongrp = $('
');\n", " var button = $('');\n", " button.click(function (evt) { fig.handle_close(fig, {}); } );\n", " button.mouseover('Stop Interaction', toolbar_mouse_event);\n", " buttongrp.append(button);\n", " var titlebar = this.root.find($('.ui-dialog-titlebar'));\n", " titlebar.prepend(buttongrp);\n", "}\n", "\n", "mpl.figure.prototype._root_extra_style = function(el){\n", " var fig = this\n", " el.on(\"remove\", function(){\n", "\tfig.close_ws(fig, {});\n", " });\n", "}\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(el){\n", " // this is important to make the div 'focusable\n", " el.attr('tabindex', 0)\n", " // reach out to IPython and tell the keyboard manager to turn it's self\n", " // off when our div gets focus\n", "\n", " // location in version 3\n", " if (IPython.notebook.keyboard_manager) {\n", " IPython.notebook.keyboard_manager.register_events(el);\n", " }\n", " else {\n", " // location in version 2\n", " IPython.keyboard_manager.register_events(el);\n", " }\n", "\n", "}\n", "\n", "mpl.figure.prototype._key_event_extra = function(event, name) {\n", " var manager = IPython.notebook.keyboard_manager;\n", " if (!manager)\n", " manager = IPython.keyboard_manager;\n", "\n", " // Check for shift+enter\n", " if (event.shiftKey && event.which == 13) {\n", " this.canvas_div.blur();\n", " event.shiftKey = false;\n", " // Send a \"J\" for go to next cell\n", " event.which = 74;\n", " event.keyCode = 74;\n", " manager.command_mode();\n", " manager.handle_keydown(event);\n", " }\n", "}\n", "\n", "mpl.figure.prototype.handle_save = function(fig, msg) {\n", " fig.ondownload(fig, null);\n", "}\n", "\n", "\n", "mpl.find_output_cell = function(html_output) {\n", " // Return the cell and output element which can be found *uniquely* in the notebook.\n", " // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n", " // IPython event is triggered only after the cells have been serialised, which for\n", " // our purposes (turning an active figure into a static one), is too late.\n", " var cells = IPython.notebook.get_cells();\n", " var ncells = cells.length;\n", " for (var i=0; i= 3 moved mimebundle to data attribute of output\n", " data = data.data;\n", " }\n", " if (data['text/html'] == html_output) {\n", " return [cell, data, j];\n", " }\n", " }\n", " }\n", " }\n", "}\n", "\n", "// Register the function which deals with the matplotlib target/channel.\n", "// The kernel may be null if the page has been refreshed.\n", "if (IPython.notebook.kernel != null) {\n", " IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n", "}\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Make the figure instance\n", "fig,ax = plt.subplots()\n", "\n", "# Plot the data\n", "ax.plot(temp,altitude[1:],'-o')\n", "\n", "# Some aestethics\n", "ax.set_xlabel('Temperature (K)')\n", "ax.set_ylabel('Altitude (m)')\n", "\n", "# To limit the vertical extent uncomment the following line\n", "#ax.set_ylim(0,3000)\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# 2D plots in time and height\n", "Now, let's combine the two previous formats. I will repeat some of the previous load statements in order to enable standalone execution of the code below." ] }, { "cell_type": "code", "execution_count": 164, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(24,)\n", "(24, 39)\n", "(24, 39)\n" ] } ], "source": [ "# Import packages\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "import pandas as pd\n", "from netCDF4 import Dataset,num2date\n", "\n", "#Don't use this command in the terminal!!!\n", "%matplotlib notebook \n", "\n", "\n", "# Make the file object\n", "wrfout = Dataset('wrfout_d03_2015-07-01_00:00:00')\n", "\n", "# Load wrf time\n", "wrftime = wrfout.variables['XTIME'][:]\n", "timeunits = wrfout.variables['XTIME'].units\n", "time = num2date(wrftime,units=timeunits,calendar='standard')\n", "\n", "# Load wrf altitude\n", "ph = wrfout.variables['PH'][:,1:,100,100] # Notice that we have now changed the vertical extent to exclude the surface\n", "phb = wrfout.variables['PHB'][:,1:,100,100] \n", "altitude = (ph+phb)/9.81\n", "\n", "# Load wrf temperature\n", "wrfT = wrfout.variables['T'][:,:,100,100]\n", "temp = wrfT+300.\n", "\n", "print time.shape\n", "print altitude.shape\n", "print temp.shape" ] }, { "cell_type": "code", "execution_count": 165, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(24, 39)\n" ] } ], "source": [ "# Make time array 2d as well\n", "twodimtime = np.tile(time, (39,1)).transpose()\n", "print twodimtime.shape" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, we have the correct shape for the date array. Note that this whole step would not have been necessary if the altitude would be constant. However, the altitude of pressure levels can change in time, and thus we will keep both dimensions. Since altitude is passed as 2d array, time needs to be 2d as well for consistency. \n", "\n", "There is one more problem with the datetime array. Matplotlib doesn't understand this format for 2d plots. Therefore, we need to do a workaround" ] }, { "cell_type": "code", "execution_count": 166, "metadata": { "collapsed": false }, "outputs": [ { "data": { "application/javascript": [ "/* Put everything inside the global mpl namespace */\n", "window.mpl = {};\n", "\n", "mpl.get_websocket_type = function() {\n", " if (typeof(WebSocket) !== 'undefined') {\n", " return WebSocket;\n", " } else if (typeof(MozWebSocket) !== 'undefined') {\n", " return MozWebSocket;\n", " } else {\n", " alert('Your browser does not have WebSocket support.' +\n", " 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n", " 'Firefox 4 and 5 are also supported but you ' +\n", " 'have to enable WebSockets in about:config.');\n", " };\n", "}\n", "\n", "mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n", " this.id = figure_id;\n", "\n", " this.ws = websocket;\n", "\n", " this.supports_binary = (this.ws.binaryType != undefined);\n", "\n", " if (!this.supports_binary) {\n", " var warnings = document.getElementById(\"mpl-warnings\");\n", " if (warnings) {\n", " warnings.style.display = 'block';\n", " warnings.textContent = (\n", " \"This browser does not support binary websocket messages. \" +\n", " \"Performance may be slow.\");\n", " }\n", " }\n", "\n", " this.imageObj = new Image();\n", "\n", " this.context = undefined;\n", " this.message = undefined;\n", " this.canvas = undefined;\n", " this.rubberband_canvas = undefined;\n", " this.rubberband_context = undefined;\n", " this.format_dropdown = undefined;\n", "\n", " this.image_mode = 'full';\n", "\n", " this.root = $('
');\n", " this._root_extra_style(this.root)\n", " this.root.attr('style', 'display: inline-block');\n", "\n", " $(parent_element).append(this.root);\n", "\n", " this._init_header(this);\n", " this._init_canvas(this);\n", " this._init_toolbar(this);\n", "\n", " var fig = this;\n", "\n", " this.waiting = false;\n", "\n", " this.ws.onopen = function () {\n", " fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n", " fig.send_message(\"send_image_mode\", {});\n", " fig.send_message(\"refresh\", {});\n", " }\n", "\n", " this.imageObj.onload = function() {\n", " if (fig.image_mode == 'full') {\n", " // Full images could contain transparency (where diff images\n", " // almost always do), so we need to clear the canvas so that\n", " // there is no ghosting.\n", " fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n", " }\n", " fig.context.drawImage(fig.imageObj, 0, 0);\n", " };\n", "\n", " this.imageObj.onunload = function() {\n", " this.ws.close();\n", " }\n", "\n", " this.ws.onmessage = this._make_on_message_function(this);\n", "\n", " this.ondownload = ondownload;\n", "}\n", "\n", "mpl.figure.prototype._init_header = function() {\n", " var titlebar = $(\n", " '
');\n", " var titletext = $(\n", " '
');\n", " titlebar.append(titletext)\n", " this.root.append(titlebar);\n", " this.header = titletext[0];\n", "}\n", "\n", "\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "\n", "mpl.figure.prototype._root_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "mpl.figure.prototype._init_canvas = function() {\n", " var fig = this;\n", "\n", " var canvas_div = $('
');\n", "\n", " canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n", "\n", " function canvas_keyboard_event(event) {\n", " return fig.key_event(event, event['data']);\n", " }\n", "\n", " canvas_div.keydown('key_press', canvas_keyboard_event);\n", " canvas_div.keyup('key_release', canvas_keyboard_event);\n", " this.canvas_div = canvas_div\n", " this._canvas_extra_style(canvas_div)\n", " this.root.append(canvas_div);\n", "\n", " var canvas = $('');\n", " canvas.addClass('mpl-canvas');\n", " canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n", "\n", " this.canvas = canvas[0];\n", " this.context = canvas[0].getContext(\"2d\");\n", "\n", " var rubberband = $('');\n", " rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n", "\n", " var pass_mouse_events = true;\n", "\n", " canvas_div.resizable({\n", " start: function(event, ui) {\n", " pass_mouse_events = false;\n", " },\n", " resize: function(event, ui) {\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " stop: function(event, ui) {\n", " pass_mouse_events = true;\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " });\n", "\n", " function mouse_event_fn(event) {\n", " if (pass_mouse_events)\n", " return fig.mouse_event(event, event['data']);\n", " }\n", "\n", " rubberband.mousedown('button_press', mouse_event_fn);\n", " rubberband.mouseup('button_release', mouse_event_fn);\n", " // Throttle sequential mouse events to 1 every 20ms.\n", " rubberband.mousemove('motion_notify', mouse_event_fn);\n", "\n", " rubberband.mouseenter('figure_enter', mouse_event_fn);\n", " rubberband.mouseleave('figure_leave', mouse_event_fn);\n", "\n", " canvas_div.on(\"wheel\", function (event) {\n", " event = event.originalEvent;\n", " event['data'] = 'scroll'\n", " if (event.deltaY < 0) {\n", " event.step = 1;\n", " } else {\n", " event.step = -1;\n", " }\n", " mouse_event_fn(event);\n", " });\n", "\n", " canvas_div.append(canvas);\n", " canvas_div.append(rubberband);\n", "\n", " this.rubberband = rubberband;\n", " this.rubberband_canvas = rubberband[0];\n", " this.rubberband_context = rubberband[0].getContext(\"2d\");\n", " this.rubberband_context.strokeStyle = \"#000000\";\n", "\n", " this._resize_canvas = function(width, height) {\n", " // Keep the size of the canvas, canvas container, and rubber band\n", " // canvas in synch.\n", " canvas_div.css('width', width)\n", " canvas_div.css('height', height)\n", "\n", " canvas.attr('width', width);\n", " canvas.attr('height', height);\n", "\n", " rubberband.attr('width', width);\n", " rubberband.attr('height', height);\n", " }\n", "\n", " // Set the figure to an initial 600x600px, this will subsequently be updated\n", " // upon first draw.\n", " this._resize_canvas(600, 600);\n", "\n", " // Disable right mouse context menu.\n", " $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n", " return false;\n", " });\n", "\n", " function set_focus () {\n", " canvas.focus();\n", " canvas_div.focus();\n", " }\n", "\n", " window.setTimeout(set_focus, 100);\n", "}\n", "\n", "mpl.figure.prototype._init_toolbar = function() {\n", " var fig = this;\n", "\n", " var nav_element = $('
')\n", " nav_element.attr('style', 'width: 100%');\n", " this.root.append(nav_element);\n", "\n", " // Define a callback function for later on.\n", " function toolbar_event(event) {\n", " return fig.toolbar_button_onclick(event['data']);\n", " }\n", " function toolbar_mouse_event(event) {\n", " return fig.toolbar_button_onmouseover(event['data']);\n", " }\n", "\n", " for(var toolbar_ind in mpl.toolbar_items) {\n", " var name = mpl.toolbar_items[toolbar_ind][0];\n", " var tooltip = mpl.toolbar_items[toolbar_ind][1];\n", " var image = mpl.toolbar_items[toolbar_ind][2];\n", " var method_name = mpl.toolbar_items[toolbar_ind][3];\n", "\n", " if (!name) {\n", " // put a spacer in here.\n", " continue;\n", " }\n", " var button = $('');\n", " button.click(method_name, toolbar_event);\n", " button.mouseover(tooltip, toolbar_mouse_event);\n", " nav_element.append(button);\n", " }\n", "\n", " // Add the status bar.\n", " var status_bar = $('');\n", " nav_element.append(status_bar);\n", " this.message = status_bar[0];\n", "\n", " // Add the close button to the window.\n", " var buttongrp = $('
');\n", " var button = $('');\n", " button.click(function (evt) { fig.handle_close(fig, {}); } );\n", " button.mouseover('Stop Interaction', toolbar_mouse_event);\n", " buttongrp.append(button);\n", " var titlebar = this.root.find($('.ui-dialog-titlebar'));\n", " titlebar.prepend(buttongrp);\n", "}\n", "\n", "mpl.figure.prototype._root_extra_style = function(el){\n", " var fig = this\n", " el.on(\"remove\", function(){\n", "\tfig.close_ws(fig, {});\n", " });\n", "}\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(el){\n", " // this is important to make the div 'focusable\n", " el.attr('tabindex', 0)\n", " // reach out to IPython and tell the keyboard manager to turn it's self\n", " // off when our div gets focus\n", "\n", " // location in version 3\n", " if (IPython.notebook.keyboard_manager) {\n", " IPython.notebook.keyboard_manager.register_events(el);\n", " }\n", " else {\n", " // location in version 2\n", " IPython.keyboard_manager.register_events(el);\n", " }\n", "\n", "}\n", "\n", "mpl.figure.prototype._key_event_extra = function(event, name) {\n", " var manager = IPython.notebook.keyboard_manager;\n", " if (!manager)\n", " manager = IPython.keyboard_manager;\n", "\n", " // Check for shift+enter\n", " if (event.shiftKey && event.which == 13) {\n", " this.canvas_div.blur();\n", " event.shiftKey = false;\n", " // Send a \"J\" for go to next cell\n", " event.which = 74;\n", " event.keyCode = 74;\n", " manager.command_mode();\n", " manager.handle_keydown(event);\n", " }\n", "}\n", "\n", "mpl.figure.prototype.handle_save = function(fig, msg) {\n", " fig.ondownload(fig, null);\n", "}\n", "\n", "\n", "mpl.find_output_cell = function(html_output) {\n", " // Return the cell and output element which can be found *uniquely* in the notebook.\n", " // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n", " // IPython event is triggered only after the cells have been serialised, which for\n", " // our purposes (turning an active figure into a static one), is too late.\n", " var cells = IPython.notebook.get_cells();\n", " var ncells = cells.length;\n", " for (var i=0; i= 3 moved mimebundle to data attribute of output\n", " data = data.data;\n", " }\n", " if (data['text/html'] == html_output) {\n", " return [cell, data, j];\n", " }\n", " }\n", " }\n", " }\n", "}\n", "\n", "// Register the function which deals with the matplotlib target/channel.\n", "// The kernel may be null if the page has been refreshed.\n", "if (IPython.notebook.kernel != null) {\n", " IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n", "}\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Create figure instance\n", "fig,ax=plt.subplots()\n", "\n", "# Create basemap instance\n", "m = Basemap(width=201*3000,height=201*3000,\\\n", " resolution='l',area_thresh=1000.,projection='lcc',\\\n", " lat_1=30.,lat_2=60.,\n", " lat_0=52.86,lon_0=3.44,ax=ax)\n", "\n", "m.drawcoastlines()\n", "\n", "# Fill background\n", "m.fillcontinents(color='coral',lake_color='aqua')\n", "#m.etopo() # requires PIL, not on hpc\n", " \n", "# Draw parallels and meridians.\n", "m.drawparallels(np.arange(-80.,81.,1.))\n", "m.drawmeridians(np.arange(-180.,181.,1.))\n", "\n", "# See what it looks like\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, we're going to project our data on the map. For this, we need to convert our longitude and latitude points (in decimal degrees) to figure coordinates (in m). " ] }, { "cell_type": "code", "execution_count": 233, "metadata": { "collapsed": false }, "outputs": [ { "data": { "application/javascript": [ "/* Put everything inside the global mpl namespace */\n", "window.mpl = {};\n", "\n", "mpl.get_websocket_type = function() {\n", " if (typeof(WebSocket) !== 'undefined') {\n", " return WebSocket;\n", " } else if (typeof(MozWebSocket) !== 'undefined') {\n", " return MozWebSocket;\n", " } else {\n", " alert('Your browser does not have WebSocket support.' +\n", " 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n", " 'Firefox 4 and 5 are also supported but you ' +\n", " 'have to enable WebSockets in about:config.');\n", " };\n", "}\n", "\n", "mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n", " this.id = figure_id;\n", "\n", " this.ws = websocket;\n", "\n", " this.supports_binary = (this.ws.binaryType != undefined);\n", "\n", " if (!this.supports_binary) {\n", " var warnings = document.getElementById(\"mpl-warnings\");\n", " if (warnings) {\n", " warnings.style.display = 'block';\n", " warnings.textContent = (\n", " \"This browser does not support binary websocket messages. \" +\n", " \"Performance may be slow.\");\n", " }\n", " }\n", "\n", " this.imageObj = new Image();\n", "\n", " this.context = undefined;\n", " this.message = undefined;\n", " this.canvas = undefined;\n", " this.rubberband_canvas = undefined;\n", " this.rubberband_context = undefined;\n", " this.format_dropdown = undefined;\n", "\n", " this.image_mode = 'full';\n", "\n", " this.root = $('
');\n", " this._root_extra_style(this.root)\n", " this.root.attr('style', 'display: inline-block');\n", "\n", " $(parent_element).append(this.root);\n", "\n", " this._init_header(this);\n", " this._init_canvas(this);\n", " this._init_toolbar(this);\n", "\n", " var fig = this;\n", "\n", " this.waiting = false;\n", "\n", " this.ws.onopen = function () {\n", " fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n", " fig.send_message(\"send_image_mode\", {});\n", " fig.send_message(\"refresh\", {});\n", " }\n", "\n", " this.imageObj.onload = function() {\n", " if (fig.image_mode == 'full') {\n", " // Full images could contain transparency (where diff images\n", " // almost always do), so we need to clear the canvas so that\n", " // there is no ghosting.\n", " fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n", " }\n", " fig.context.drawImage(fig.imageObj, 0, 0);\n", " };\n", "\n", " this.imageObj.onunload = function() {\n", " this.ws.close();\n", " }\n", "\n", " this.ws.onmessage = this._make_on_message_function(this);\n", "\n", " this.ondownload = ondownload;\n", "}\n", "\n", "mpl.figure.prototype._init_header = function() {\n", " var titlebar = $(\n", " '
');\n", " var titletext = $(\n", " '
');\n", " titlebar.append(titletext)\n", " this.root.append(titlebar);\n", " this.header = titletext[0];\n", "}\n", "\n", "\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "\n", "mpl.figure.prototype._root_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "mpl.figure.prototype._init_canvas = function() {\n", " var fig = this;\n", "\n", " var canvas_div = $('
');\n", "\n", " canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n", "\n", " function canvas_keyboard_event(event) {\n", " return fig.key_event(event, event['data']);\n", " }\n", "\n", " canvas_div.keydown('key_press', canvas_keyboard_event);\n", " canvas_div.keyup('key_release', canvas_keyboard_event);\n", " this.canvas_div = canvas_div\n", " this._canvas_extra_style(canvas_div)\n", " this.root.append(canvas_div);\n", "\n", " var canvas = $('');\n", " canvas.addClass('mpl-canvas');\n", " canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n", "\n", " this.canvas = canvas[0];\n", " this.context = canvas[0].getContext(\"2d\");\n", "\n", " var rubberband = $('');\n", " rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n", "\n", " var pass_mouse_events = true;\n", "\n", " canvas_div.resizable({\n", " start: function(event, ui) {\n", " pass_mouse_events = false;\n", " },\n", " resize: function(event, ui) {\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " stop: function(event, ui) {\n", " pass_mouse_events = true;\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " });\n", "\n", " function mouse_event_fn(event) {\n", " if (pass_mouse_events)\n", " return fig.mouse_event(event, event['data']);\n", " }\n", "\n", " rubberband.mousedown('button_press', mouse_event_fn);\n", " rubberband.mouseup('button_release', mouse_event_fn);\n", " // Throttle sequential mouse events to 1 every 20ms.\n", " rubberband.mousemove('motion_notify', mouse_event_fn);\n", "\n", " rubberband.mouseenter('figure_enter', mouse_event_fn);\n", " rubberband.mouseleave('figure_leave', mouse_event_fn);\n", "\n", " canvas_div.on(\"wheel\", function (event) {\n", " event = event.originalEvent;\n", " event['data'] = 'scroll'\n", " if (event.deltaY < 0) {\n", " event.step = 1;\n", " } else {\n", " event.step = -1;\n", " }\n", " mouse_event_fn(event);\n", " });\n", "\n", " canvas_div.append(canvas);\n", " canvas_div.append(rubberband);\n", "\n", " this.rubberband = rubberband;\n", " this.rubberband_canvas = rubberband[0];\n", " this.rubberband_context = rubberband[0].getContext(\"2d\");\n", " this.rubberband_context.strokeStyle = \"#000000\";\n", "\n", " this._resize_canvas = function(width, height) {\n", " // Keep the size of the canvas, canvas container, and rubber band\n", " // canvas in synch.\n", " canvas_div.css('width', width)\n", " canvas_div.css('height', height)\n", "\n", " canvas.attr('width', width);\n", " canvas.attr('height', height);\n", "\n", " rubberband.attr('width', width);\n", " rubberband.attr('height', height);\n", " }\n", "\n", " // Set the figure to an initial 600x600px, this will subsequently be updated\n", " // upon first draw.\n", " this._resize_canvas(600, 600);\n", "\n", " // Disable right mouse context menu.\n", " $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n", " return false;\n", " });\n", "\n", " function set_focus () {\n", " canvas.focus();\n", " canvas_div.focus();\n", " }\n", "\n", " window.setTimeout(set_focus, 100);\n", "}\n", "\n", "mpl.figure.prototype._init_toolbar = function() {\n", " var fig = this;\n", "\n", " var nav_element = $('
')\n", " nav_element.attr('style', 'width: 100%');\n", " this.root.append(nav_element);\n", "\n", " // Define a callback function for later on.\n", " function toolbar_event(event) {\n", " return fig.toolbar_button_onclick(event['data']);\n", " }\n", " function toolbar_mouse_event(event) {\n", " return fig.toolbar_button_onmouseover(event['data']);\n", " }\n", "\n", " for(var toolbar_ind in mpl.toolbar_items) {\n", " var name = mpl.toolbar_items[toolbar_ind][0];\n", " var tooltip = mpl.toolbar_items[toolbar_ind][1];\n", " var image = mpl.toolbar_items[toolbar_ind][2];\n", " var method_name = mpl.toolbar_items[toolbar_ind][3];\n", "\n", " if (!name) {\n", " // put a spacer in here.\n", " continue;\n", " }\n", " var button = $('');\n", " button.click(method_name, toolbar_event);\n", " button.mouseover(tooltip, toolbar_mouse_event);\n", " nav_element.append(button);\n", " }\n", "\n", " // Add the status bar.\n", " var status_bar = $('');\n", " nav_element.append(status_bar);\n", " this.message = status_bar[0];\n", "\n", " // Add the close button to the window.\n", " var buttongrp = $('
');\n", " var button = $('');\n", " button.click(function (evt) { fig.handle_close(fig, {}); } );\n", " button.mouseover('Stop Interaction', toolbar_mouse_event);\n", " buttongrp.append(button);\n", " var titlebar = this.root.find($('.ui-dialog-titlebar'));\n", " titlebar.prepend(buttongrp);\n", "}\n", "\n", "mpl.figure.prototype._root_extra_style = function(el){\n", " var fig = this\n", " el.on(\"remove\", function(){\n", "\tfig.close_ws(fig, {});\n", " });\n", "}\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(el){\n", " // this is important to make the div 'focusable\n", " el.attr('tabindex', 0)\n", " // reach out to IPython and tell the keyboard manager to turn it's self\n", " // off when our div gets focus\n", "\n", " // location in version 3\n", " if (IPython.notebook.keyboard_manager) {\n", " IPython.notebook.keyboard_manager.register_events(el);\n", " }\n", " else {\n", " // location in version 2\n", " IPython.keyboard_manager.register_events(el);\n", " }\n", "\n", "}\n", "\n", "mpl.figure.prototype._key_event_extra = function(event, name) {\n", " var manager = IPython.notebook.keyboard_manager;\n", " if (!manager)\n", " manager = IPython.keyboard_manager;\n", "\n", " // Check for shift+enter\n", " if (event.shiftKey && event.which == 13) {\n", " this.canvas_div.blur();\n", " event.shiftKey = false;\n", " // Send a \"J\" for go to next cell\n", " event.which = 74;\n", " event.keyCode = 74;\n", " manager.command_mode();\n", " manager.handle_keydown(event);\n", " }\n", "}\n", "\n", "mpl.figure.prototype.handle_save = function(fig, msg) {\n", " fig.ondownload(fig, null);\n", "}\n", "\n", "\n", "mpl.find_output_cell = function(html_output) {\n", " // Return the cell and output element which can be found *uniquely* in the notebook.\n", " // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n", " // IPython event is triggered only after the cells have been serialised, which for\n", " // our purposes (turning an active figure into a static one), is too late.\n", " var cells = IPython.notebook.get_cells();\n", " var ncells = cells.length;\n", " for (var i=0; i= 3 moved mimebundle to data attribute of output\n", " data = data.data;\n", " }\n", " if (data['text/html'] == html_output) {\n", " return [cell, data, j];\n", " }\n", " }\n", " }\n", " }\n", "}\n", "\n", "// Register the function which deals with the matplotlib target/channel.\n", "// The kernel may be null if the page has been refreshed.\n", "if (IPython.notebook.kernel != null) {\n", " IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n", "}\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Import packages\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "from netCDF4 import Dataset,num2date\n", "from mpl_toolkits.basemap import Basemap\n", "#Don't use this command in the terminal!!!\n", "%matplotlib notebook \n", "\n", "timestep = 12\n", "\n", "# Make the file object\n", "wrfin = Dataset('wrfinput_d03')\n", "wrfout = Dataset('wrfout_d03_2015-07-01_00:00:00')\n", "\n", "wrflongitude = wrfin.variables['XLONG'][:]\n", "wrflatitude = wrfin.variables['XLAT'][:]\n", "\n", "# To get rid of the extra dimension in the grid (which contains only 1 record), use numpy.squeeze\n", "longitude = wrflongitude.squeeze()\n", "latitude = wrflatitude.squeeze()\n", "\n", "# Load wrf 2m temperature (this is a seperate variable) and 10 m wind\n", "# Use fixed time and surface level\n", "temp2m = wrfout.variables['T2'][timestep,:,:]\n", "u10 = wrfout.variables['U10'][timestep,:,:]\n", "v10 = wrfout.variables['V10'][timestep,:,:]\n", "\n", "# Lets put the basemap on a figure instance\n", "fig, ax = plt.subplots()\n", "\n", "# Let's use mean lat and lon as reference point\n", "reflat = latitude.mean()\n", "reflon = longitude.mean()\n", "\n", "# Create basemap instance again, because we already used plt.show()\n", "m = Basemap(width=201*3000,height=201*3000,\\\n", " resolution='l',area_thresh=1000.,projection='lcc',\\\n", " lat_1=30.,lat_2=60.,\n", " lat_0=reflat,lon_0=reflon,ax=ax)\n", "\n", "# Convert lat/lon to map coordinates\n", "x,y = m(longitude,latitude)\n", "\n", "# Plot the data for 2m temperature\n", "t2mplot = m.pcolormesh(x,y,temp2m)\n", "\n", "# For the quiver arrows, we don't want every singly point to be plotted, because then we get 201x201 arrows.\n", "# Therefore, use a smart resample method:\n", "skip=(slice(None,None,10),slice(None,None,10))\n", "m.quiver(x[skip],y[skip],u10[skip],v10[skip])\n", "# Now, every 10th point in both directions is plotted.\n", "\n", "# Some orientation lines\n", "m.drawcoastlines()\n", "m.drawparallels(np.arange(-80.,81.,2.)) \n", "m.drawmeridians(np.arange(-180.,181.,2.))\n", "\n", "plt.colorbar(t2mplot)\n", "# Show the figure\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Note**: for 10m wind, the wrfoutput has been automatically regridded to center points. For other levels of wind, they are defined at the cell boundaries. Therefore, U and V have different shapes, and some more advanced resampling is needed before plotting, as shown below" ] }, { "cell_type": "code", "execution_count": 237, "metadata": { "collapsed": true }, "outputs": [], "source": [ "# Need to resample the u and v values so that they coincide on the grid\n", "#\n", "# --v-- --v-- \n", "# | | |\n", "# u q u q u\n", "# | | |\n", "# --v-- --v--\n", "# | | |\n", "# u q u q u\n", "# | | |\n", "# --v-- --v--\n", "#\n", "# We want to evaluate v and u simulatenously at the latitude and longitude of q\n", "# u is oversampled in the x-direction (axis = 1)\n", "# v is oversampled in the y-direction (axis = 0)\n", "\n", "# These functions should do the trick, but only for 2d arrays!!\n", "def resampleU(inputarray):\n", " i,j = inputarray.shape #works only for 2d arrays!\n", " output = np.zeros((i,j-1))\n", " for a in range(j-1):\n", " output[:,a] = (inputarray[:,a]+inputarray[:,a+1])/2.\n", " return output\n", "\n", "def resampleV(inputarray):\n", " i,j = inputarray.shape #works only for 2d arrays!\n", " output = np.zeros((i-1,j))\n", " for a in range(i-1):\n", " output[a,:] = (inputarray[a,:]+inputarray[a+1,:])/2.\n", " return output\n", "\n", "u_center = resampleU(uwrf)\n", "v_center = resampleV(vwrf)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# How to find your gridpoint\n", "In the above examples for timeseries and vertical profiles, I've just taken a random grid point. Let's say we want the gridpoint belonging to the coordinates of Cabauw: 51.971N 4.925E. Then we could use the following code:" ] }, { "cell_type": "code", "execution_count": 324, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Index: (80, 146)\n", "Coordinates: (4.9724426, 51.963196)\n" ] } ], "source": [ "def findpoint(gridlons,gridlats,point):\n", " ''' Return the index of the grid point nearest to a given point\n", " gridlons and gridlats should be 2d meshgrid arrays\n", " point is a tuple with (lon,lat)'''\n", " # First, convert long and lat arrays to flat arrays with ravel\n", " # Then, zip the two to get an array with tuples for each point on the grid\n", " grid = zip(gridlons.ravel(),gridlats.ravel())\n", " \n", " # The euclidean distance can be computed with linalg.norm \n", " euclidean = lambda a,b: np.linalg.norm(np.array(a)-np.array(b))\n", " \n", " # Compute the euclidian distance between each grid point and the query point\n", " dist = np.asarray([euclidean(cell,point) for cell in grid])\n", " \n", " # Find the minimum of the euclidian distances array and return its index\n", " # Use np.unravel_index to transform back to 2d index instead of the flattened array\n", " ind = np.unravel_index(dist.argmin(),gridlons.shape)\n", " \n", " return ind,gridlons[ind],gridlats[ind]\n", "\n", "index,ilat,ilon = findpoint(longitude,latitude,(4.971,51.971))\n", "print 'Index:',index\n", "print 'Coordinates:',(ilat,ilon)\n", " " ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 2", "language": "python", "name": "python2" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 2 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython2", "version": "2.7.11" } }, "nbformat": 4, "nbformat_minor": 0 }