<!DOCTYPE html>
<html class="writer-html5" lang="en" data-content_root="../../../">
<head>
<meta charset="utf-8" />
<meta name="viewport" content="width=device-width, initial-scale=1.0" />
<title>opendrift.readers.reader_global_landmask — OpenDrift documentation</title>
<link rel="stylesheet" type="text/css" href="../../../_static/pygments.css?v=b86133f3" />
<link rel="stylesheet" type="text/css" href="../../../_static/css/theme.css?v=e59714d7" />
<link rel="stylesheet" type="text/css" href="../../../_static/graphviz.css?v=4ae1632d" />
<link rel="stylesheet" type="text/css" href="../../../_static/plot_directive.css" />
<link rel="stylesheet" type="text/css" href="../../../_static/sg_gallery.css?v=d2d258e8" />
<link rel="stylesheet" type="text/css" href="../../../_static/sg_gallery-binder.css?v=f4aeca0c" />
<link rel="stylesheet" type="text/css" href="../../../_static/sg_gallery-dataframe.css?v=2082cf3c" />
<link rel="stylesheet" type="text/css" href="../../../_static/sg_gallery-rendered-html.css?v=1277b6f3" />
<link rel="stylesheet" type="text/css" href="../../../_static/theme_overrides.css" />
<script src="../../../_static/jquery.js?v=5d32c60e"></script>
<script src="../../../_static/_sphinx_javascript_frameworks_compat.js?v=2cd50e6c"></script>
<script src="../../../_static/documentation_options.js?v=5929fcd5"></script>
<script src="../../../_static/doctools.js?v=9bcbadda"></script>
<script src="../../../_static/sphinx_highlight.js?v=dc90522c"></script>
<script src="../../../_static/js/theme.js"></script>
<link rel="index" title="Index" href="../../../genindex.html" />
<link rel="search" title="Search" href="../../../search.html" />
</head>
<body class="wy-body-for-nav">
<div class="wy-grid-for-nav">
<nav data-toggle="wy-nav-shift" class="wy-nav-side">
<div class="wy-side-scroll">
<div class="wy-side-nav-search" >
<a href="../../../index.html" class="icon icon-home">
OpenDrift
<img src="../../../_static/opendrift_logo.png" class="logo" alt="Logo"/>
</a>
<div role="search">
<form id="rtd-search-form" class="wy-form" action="../../../search.html" method="get">
<input type="text" name="q" placeholder="Search docs" aria-label="Search docs" />
<input type="hidden" name="check_keywords" value="yes" />
<input type="hidden" name="area" value="default" />
</form>
</div>
</div><div class="wy-menu wy-menu-vertical" data-spy="affix" role="navigation" aria-label="Navigation menu">
<ul>
<li class="toctree-l1"><a class="reference internal" href="../../../index.html">Introduction to OpenDrift</a></li>
</ul>
<ul>
<li class="toctree-l1"><a class="reference internal" href="../../../history_link.html">History</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../install.html">Installing OpenDrift</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../performance.html">Performance in OpenDrift</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../tutorial.html">Tutorial</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../theory/index.html">Theory</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../theory/index.html#drift-in-the-ocean">Drift in the Ocean</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../choosing_a_model.html">How to choose which model to use</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../writing_a_new_model.html">How to write a new module</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../gallery/index.html">Gallery</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../oil_types.html">Oil types</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../interaction_with_coastline.html">Interaction with coastline</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../docker.html">Using OpenDrift in a container</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../gui.html">Graphical User Interface</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../references.html">Publications</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../services.html">Services using OpenDrift</a></li>
<li class="toctree-l1"><a class="reference internal" href="../../../autoapi/index.html">API Reference</a></li>
</ul>
</div>
</div>
</nav>
<section data-toggle="wy-nav-shift" class="wy-nav-content-wrap"><nav class="wy-nav-top" aria-label="Mobile navigation menu" >
<i data-toggle="wy-nav-top" class="fa fa-bars"></i>
<a href="../../../index.html">OpenDrift</a>
</nav>
<div class="wy-nav-content">
<div class="rst-content">
<div role="navigation" aria-label="Page navigation">
<ul class="wy-breadcrumbs">
<li><a href="../../../index.html" class="icon icon-home" aria-label="Home"></a></li>
<li class="breadcrumb-item"><a href="../../index.html">Module code</a></li>
<li class="breadcrumb-item"><a href="../../opendrift.html">opendrift</a></li>
<li class="breadcrumb-item"><a href="../readers.html">opendrift.readers</a></li>
<li class="breadcrumb-item active">opendrift.readers.reader_global_landmask</li>
<li class="wy-breadcrumbs-aside">
</li>
</ul>
<hr/>
</div>
<div role="main" class="document" itemscope="itemscope" itemtype="http://schema.org/Article">
<div itemprop="articleBody">
<h1>Source code for opendrift.readers.reader_global_landmask</h1><div class="highlight"><pre>
<span></span><span class="c1"># This file is part of OpenDrift.</span>
<span class="c1">#</span>
<span class="c1"># OpenDrift is free software: you can redistribute it and/or modify</span>
<span class="c1"># it under the terms of the GNU General Public License as published by</span>
<span class="c1"># the Free Software Foundation, version 2</span>
<span class="c1">#</span>
<span class="c1"># OpenDrift is distributed in the hope that it will be useful,</span>
<span class="c1"># but WITHOUT ANY WARRANTY; without even the implied warranty of</span>
<span class="c1"># MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the</span>
<span class="c1"># GNU General Public License for more details.</span>
<span class="c1">#</span>
<span class="c1"># You should have received a copy of the GNU General Public License</span>
<span class="c1"># along with OpenDrift. If not, see <https://www.gnu.org/licenses/>.</span>
<span class="c1">#</span>
<span class="c1"># Copyright 2020, Gaute Hope, MET Norway</span>
<span class="kn">from</span><span class="w"> </span><span class="nn">opendrift.readers.basereader</span><span class="w"> </span><span class="kn">import</span> <span class="n">BaseReader</span><span class="p">,</span> <span class="n">ContinuousReader</span>
<span class="kn">import</span><span class="w"> </span><span class="nn">pyproj</span>
<span class="kn">import</span><span class="w"> </span><span class="nn">numpy</span><span class="w"> </span><span class="k">as</span><span class="w"> </span><span class="nn">np</span>
<span class="kn">import</span><span class="w"> </span><span class="nn">cartopy.crs</span><span class="w"> </span><span class="k">as</span><span class="w"> </span><span class="nn">ccrs</span>
<span class="kn">import</span><span class="w"> </span><span class="nn">cartopy.feature</span><span class="w"> </span><span class="k">as</span><span class="w"> </span><span class="nn">cfeature</span>
<span class="kn">import</span><span class="w"> </span><span class="nn">cartopy.io.shapereader</span><span class="w"> </span><span class="k">as</span><span class="w"> </span><span class="nn">shapereader</span>
<span class="kn">import</span><span class="w"> </span><span class="nn">shapely.geometry</span><span class="w"> </span><span class="k">as</span><span class="w"> </span><span class="nn">sgeom</span>
<span class="kn">from</span><span class="w"> </span><span class="nn">shapely</span><span class="w"> </span><span class="kn">import</span> <span class="n">wkb</span><span class="p">,</span> <span class="n">clip_by_rect</span>
<span class="kn">import</span><span class="w"> </span><span class="nn">logging</span>
<span class="kn">import</span><span class="w"> </span><span class="nn">roaring_landmask</span>
<span class="n">logger</span> <span class="o">=</span> <span class="n">logging</span><span class="o">.</span><span class="n">getLogger</span><span class="p">(</span><span class="vm">__name__</span><span class="p">)</span>
<span class="c1"># TODO: replace this with weakref + thread-safety</span>
<span class="n">__roaring_mask__</span> <span class="o">=</span> <span class="kc">None</span>
<span class="n">__polys__</span> <span class="o">=</span> <span class="kc">None</span>
<div class="viewcode-block" id="get_mask">
<a class="viewcode-back" href="../../../autoapi/opendrift/readers/reader_global_landmask/index.html#opendrift.readers.reader_global_landmask.get_mask">[docs]</a>
<span class="k">def</span><span class="w"> </span><span class="nf">get_mask</span><span class="p">():</span>
<span class="w"> </span><span class="sd">"""</span>
<span class="sd"> Returns an instance of the landmask type and landmask. The mask data is</span>
<span class="sd"> usually shared between threads.</span>
<span class="sd"> """</span>
<span class="k">global</span> <span class="n">__roaring_mask__</span>
<span class="k">if</span> <span class="n">__roaring_mask__</span> <span class="ow">is</span> <span class="kc">None</span><span class="p">:</span>
<span class="kn">from</span><span class="w"> </span><span class="nn">roaring_landmask</span><span class="w"> </span><span class="kn">import</span> <span class="n">RoaringLandmask</span>
<span class="n">__roaring_mask__</span> <span class="o">=</span> <span class="n">RoaringLandmask</span><span class="o">.</span><span class="n">new</span><span class="p">()</span>
<span class="k">return</span> <span class="n">__roaring_mask__</span></div>
<div class="viewcode-block" id="LandmaskFeature">
<a class="viewcode-back" href="../../../autoapi/opendrift/readers/reader_global_landmask/index.html#opendrift.readers.reader_global_landmask.LandmaskFeature">[docs]</a>
<span class="k">class</span><span class="w"> </span><span class="nc">LandmaskFeature</span><span class="p">(</span><span class="n">cfeature</span><span class="o">.</span><span class="n">GSHHSFeature</span><span class="p">):</span>
<span class="k">def</span><span class="w"> </span><span class="fm">__init__</span><span class="p">(</span><span class="bp">self</span><span class="p">,</span> <span class="n">scale</span><span class="o">=</span><span class="s1">'auto'</span><span class="p">,</span> <span class="n">globe</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span> <span class="o">**</span><span class="n">kwargs</span><span class="p">):</span>
<span class="nb">super</span><span class="p">()</span><span class="o">.</span><span class="fm">__init__</span><span class="p">(</span><span class="n">scale</span><span class="p">,</span> <span class="o">**</span><span class="n">kwargs</span><span class="p">)</span>
<span class="k">if</span> <span class="n">globe</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span><span class="p">:</span>
<span class="bp">self</span><span class="o">.</span><span class="n">_crs</span> <span class="o">=</span> <span class="n">ccrs</span><span class="o">.</span><span class="n">PlateCarree</span><span class="p">(</span><span class="n">globe</span><span class="o">=</span><span class="n">globe</span><span class="p">)</span>
<div class="viewcode-block" id="LandmaskFeature.geometries">
<a class="viewcode-back" href="../../../autoapi/opendrift/readers/reader_global_landmask/index.html#opendrift.readers.reader_global_landmask.LandmaskFeature.geometries">[docs]</a>
<span class="k">def</span><span class="w"> </span><span class="nf">geometries</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
<span class="bp">self</span><span class="o">.</span><span class="n">intersecting_geometries</span><span class="p">(</span><span class="n">extent</span><span class="o">=</span><span class="kc">None</span><span class="p">)</span></div>
<div class="viewcode-block" id="LandmaskFeature.intersecting_geometries">
<a class="viewcode-back" href="../../../autoapi/opendrift/readers/reader_global_landmask/index.html#opendrift.readers.reader_global_landmask.LandmaskFeature.intersecting_geometries">[docs]</a>
<span class="k">def</span><span class="w"> </span><span class="nf">intersecting_geometries</span><span class="p">(</span><span class="bp">self</span><span class="p">,</span> <span class="n">extent</span><span class="p">):</span>
<span class="k">if</span> <span class="bp">self</span><span class="o">.</span><span class="n">_scale</span> <span class="o">==</span> <span class="s1">'auto'</span><span class="p">:</span>
<span class="n">scale</span> <span class="o">=</span> <span class="bp">self</span><span class="o">.</span><span class="n">_scale_from_extent</span><span class="p">(</span><span class="n">extent</span><span class="p">)</span>
<span class="k">else</span><span class="p">:</span>
<span class="n">scale</span> <span class="o">=</span> <span class="bp">self</span><span class="o">.</span><span class="n">_scale</span><span class="p">[</span><span class="mi">0</span><span class="p">]</span>
<span class="k">if</span> <span class="n">extent</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span><span class="p">:</span>
<span class="n">extent_geom</span> <span class="o">=</span> <span class="n">sgeom</span><span class="o">.</span><span class="n">box</span><span class="p">(</span><span class="n">extent</span><span class="p">[</span><span class="mi">0</span><span class="p">],</span> <span class="n">extent</span><span class="p">[</span><span class="mi">2</span><span class="p">],</span>
<span class="n">extent</span><span class="p">[</span><span class="mi">1</span><span class="p">],</span> <span class="n">extent</span><span class="p">[</span><span class="mi">3</span><span class="p">])</span>
<span class="k">for</span> <span class="n">level</span> <span class="ow">in</span> <span class="bp">self</span><span class="o">.</span><span class="n">_levels</span><span class="p">:</span>
<span class="n">geoms</span> <span class="o">=</span> <span class="n">cfeature</span><span class="o">.</span><span class="n">GSHHSFeature</span><span class="o">.</span><span class="n">_geometries_cache</span><span class="o">.</span><span class="n">get</span><span class="p">((</span><span class="n">scale</span><span class="p">,</span> <span class="n">level</span><span class="p">,</span> <span class="n">extent_geom</span><span class="p">))</span>
<span class="k">if</span> <span class="n">geoms</span> <span class="ow">is</span> <span class="kc">None</span><span class="p">:</span>
<span class="k">if</span> <span class="n">scale</span> <span class="ow">in</span> <span class="p">[</span><span class="s1">'f'</span><span class="p">,</span> <span class="s1">'full'</span><span class="p">]:</span> <span class="c1"># Getting fullres polygons from roaring_landmask</span>
<span class="n">geoms</span> <span class="o">=</span> <span class="n">cfeature</span><span class="o">.</span><span class="n">GSHHSFeature</span><span class="o">.</span><span class="n">_geometries_cache</span><span class="o">.</span><span class="n">get</span><span class="p">(</span><span class="s1">'ROARING'</span><span class="p">)</span>
<span class="k">if</span> <span class="n">geoms</span> <span class="ow">is</span> <span class="kc">None</span><span class="p">:</span>
<span class="n">logger</span><span class="o">.</span><span class="n">debug</span><span class="p">(</span><span class="s1">'Getting fullres shapes from roaring landmask'</span><span class="p">)</span>
<span class="n">provider</span> <span class="o">=</span> <span class="n">roaring_landmask</span><span class="o">.</span><span class="n">LandmaskProvider</span><span class="o">.</span><span class="n">Gshhg</span>
<span class="n">shapes</span> <span class="o">=</span> <span class="n">roaring_landmask</span><span class="o">.</span><span class="n">Shapes</span><span class="o">.</span><span class="n">new</span><span class="p">(</span><span class="n">provider</span><span class="p">)</span>
<span class="n">polys</span> <span class="o">=</span> <span class="n">roaring_landmask</span><span class="o">.</span><span class="n">Shapes</span><span class="o">.</span><span class="n">wkb</span><span class="p">(</span><span class="n">provider</span><span class="p">)</span>
<span class="n">__polys__</span> <span class="o">=</span> <span class="n">wkb</span><span class="o">.</span><span class="n">loads</span><span class="p">(</span><span class="n">polys</span><span class="p">)</span>
<span class="n">geoms</span> <span class="o">=</span> <span class="n">__polys__</span><span class="o">.</span><span class="n">geoms</span>
<span class="n">cfeature</span><span class="o">.</span><span class="n">GSHHSFeature</span><span class="o">.</span><span class="n">_geometries_cache</span><span class="p">[</span><span class="s1">'ROARING'</span><span class="p">]</span> <span class="o">=</span> <span class="n">geoms</span>
<span class="k">else</span><span class="p">:</span>
<span class="n">logger</span><span class="o">.</span><span class="n">debug</span><span class="p">(</span><span class="s1">'Getting fullres shapes from roaring landmask cache'</span><span class="p">)</span>
<span class="k">else</span><span class="p">:</span>
<span class="n">logger</span><span class="o">.</span><span class="n">debug</span><span class="p">(</span><span class="sa">f</span><span class="s1">'Loading shapes (</span><span class="se">\'</span><span class="si">{</span><span class="n">scale</span><span class="si">}</span><span class="se">\'</span><span class="s1"> level </span><span class="si">{</span><span class="n">level</span><span class="si">}</span><span class="s1">) with Cartopy shapereader...'</span><span class="p">)</span>
<span class="c1"># Load GSHHS geometries from appropriate shape file.</span>
<span class="c1"># TODO selective load based on bbox of each geom in file.</span>
<span class="n">path</span> <span class="o">=</span> <span class="n">shapereader</span><span class="o">.</span><span class="n">gshhs</span><span class="p">(</span><span class="n">scale</span><span class="p">,</span> <span class="n">level</span><span class="p">)</span>
<span class="n">geoms</span> <span class="o">=</span> <span class="nb">tuple</span><span class="p">(</span><span class="n">shapereader</span><span class="o">.</span><span class="n">Reader</span><span class="p">(</span><span class="n">path</span><span class="p">)</span><span class="o">.</span><span class="n">geometries</span><span class="p">())</span>
<span class="c1"># Unlike Caropy.GSHHSFeature, we clip the polygons to plot extent</span>
<span class="c1">#geoms = [g.intersection(extent_geom) for g in geoms if g.intersects(extent_geom)]</span>
<span class="n">delta</span> <span class="o">=</span> <span class="mf">.1</span> <span class="c1"># Adding small delta to avoid segments along map borders</span>
<span class="n">geoms</span> <span class="o">=</span> <span class="p">[</span><span class="n">clip_by_rect</span><span class="p">(</span><span class="n">g</span><span class="p">,</span> <span class="n">extent</span><span class="p">[</span><span class="mi">0</span><span class="p">]</span><span class="o">-</span><span class="n">delta</span><span class="p">,</span> <span class="n">extent</span><span class="p">[</span><span class="mi">2</span><span class="p">]</span><span class="o">-</span><span class="n">delta</span><span class="p">,</span> <span class="n">extent</span><span class="p">[</span><span class="mi">1</span><span class="p">]</span><span class="o">+</span><span class="n">delta</span><span class="p">,</span> <span class="n">extent</span><span class="p">[</span><span class="mi">3</span><span class="p">]</span><span class="o">+</span><span class="n">delta</span><span class="p">)</span>
<span class="k">for</span> <span class="n">g</span> <span class="ow">in</span> <span class="n">geoms</span> <span class="k">if</span> <span class="n">g</span><span class="o">.</span><span class="n">intersects</span><span class="p">(</span><span class="n">extent_geom</span><span class="p">)]</span>
<span class="c1"># Due to the above, the cache is also depending on extent</span>
<span class="c1"># Note: if plotting several extents within same session, it would be benefitial</span>
<span class="c1"># with another cache independent of extent</span>
<span class="n">cfeature</span><span class="o">.</span><span class="n">GSHHSFeature</span><span class="o">.</span><span class="n">_geometries_cache</span><span class="p">[(</span><span class="n">scale</span><span class="p">,</span> <span class="n">level</span><span class="p">,</span> <span class="n">extent_geom</span><span class="p">)]</span> <span class="o">=</span> <span class="n">geoms</span>
<span class="k">for</span> <span class="n">geom</span> <span class="ow">in</span> <span class="n">geoms</span><span class="p">:</span>
<span class="k">yield</span> <span class="n">geom</span>
<span class="k">if</span> <span class="n">scale</span> <span class="ow">in</span> <span class="p">[</span><span class="s1">'f'</span><span class="p">,</span> <span class="s1">'full'</span><span class="p">]:</span>
<span class="k">return</span> <span class="c1"># Only a single (combined) level with full resolution troaring_landmask</span>
<span class="c1"># # If scale is full use the geometries from roaring landmask, otherwise</span>
<span class="c1"># # fall back to Cartopy provider.</span>
<span class="c1"># if scale == 'f':</span>
<span class="c1"># logger.debug(f"Adding full GSHHG shapes from roaring-landmask, extent: {extent}..")</span>
<span class="c1"># if __polys__ is None:</span>
<span class="c1"># from roaring_landmask import Gshhg</span>
<span class="c1"># polys = Gshhg.wkb()</span>
<span class="c1"># __polys__ = wkb.loads(polys)</span>
<span class="c1"># if extent is not None:</span>
<span class="c1"># extent = box(extent[0], extent[2], extent[1], extent[3])</span>
<span class="c1"># extent = shapely.prepared.prep(extent)</span>
<span class="c1"># return (p for p in __polys__.geoms if extent.intersects(p))</span>
<span class="c1"># else:</span>
<span class="c1"># return __polys__.geoms</span>
<span class="c1"># else:</span>
<span class="n">logger</span><span class="o">.</span><span class="n">debug</span><span class="p">(</span><span class="sa">f</span><span class="s2">"Adding GSHHG shapes from cartopy, scale: </span><span class="si">{</span><span class="n">scale</span><span class="si">}</span><span class="s2">, extent: </span><span class="si">{</span><span class="n">extent</span><span class="si">}</span><span class="s2">.."</span><span class="p">)</span>
<span class="k">return</span> <span class="nb">super</span><span class="p">()</span><span class="o">.</span><span class="n">intersecting_geometries</span><span class="p">(</span><span class="n">extent</span><span class="p">)</span></div>
</div>
<div class="viewcode-block" id="plot_land">
<a class="viewcode-back" href="../../../autoapi/opendrift/readers/reader_global_landmask/index.html#opendrift.readers.reader_global_landmask.plot_land">[docs]</a>
<span class="k">def</span><span class="w"> </span><span class="nf">plot_land</span><span class="p">(</span><span class="n">ax</span><span class="p">,</span> <span class="n">lonmin</span><span class="p">,</span> <span class="n">latmin</span><span class="p">,</span> <span class="n">lonmax</span><span class="p">,</span> <span class="n">latmax</span><span class="p">,</span> <span class="n">fast</span><span class="p">,</span> <span class="n">ocean_color</span> <span class="o">=</span> <span class="s1">'white'</span><span class="p">,</span> <span class="n">land_color</span> <span class="o">=</span> <span class="n">cfeature</span><span class="o">.</span><span class="n">COLORS</span><span class="p">[</span><span class="s1">'land'</span><span class="p">],</span> <span class="n">lscale</span> <span class="o">=</span> <span class="s1">'auto'</span><span class="p">,</span> <span class="n">crs_plot</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span> <span class="n">crs_lonlat</span><span class="o">=</span><span class="kc">None</span><span class="p">):</span>
<span class="w"> </span><span class="sd">"""</span>
<span class="sd"> Plot the landmask or the shapes from GSHHG.</span>
<span class="sd"> """</span>
<span class="k">def</span><span class="w"> </span><span class="nf">show_landmask_roaring</span><span class="p">(</span><span class="n">roaring</span><span class="p">):</span>
<span class="n">maxn</span> <span class="o">=</span> <span class="mf">512.</span>
<span class="k">if</span> <span class="n">crs_plot</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span> <span class="ow">and</span> <span class="n">crs_lonlat</span> <span class="ow">is</span> <span class="ow">not</span> <span class="kc">None</span><span class="p">:</span>
<span class="c1"># Make transformer and convert lon,lat to Mercator x,y</span>
<span class="n">transformer</span> <span class="o">=</span> <span class="n">pyproj</span><span class="o">.</span><span class="n">Transformer</span><span class="o">.</span><span class="n">from_crs</span><span class="p">(</span><span class="n">crs_lonlat</span><span class="p">,</span> <span class="n">crs_plot</span><span class="p">)</span>
<span class="n">xmin</span><span class="p">,</span> <span class="n">ymin</span> <span class="o">=</span> <span class="n">transformer</span><span class="o">.</span><span class="n">transform</span><span class="p">(</span><span class="n">lonmin</span><span class="p">,</span> <span class="n">latmin</span><span class="p">)</span>
<span class="n">xmax</span><span class="p">,</span> <span class="n">ymax</span> <span class="o">=</span> <span class="n">transformer</span><span class="o">.</span><span class="n">transform</span><span class="p">(</span><span class="n">lonmax</span><span class="p">,</span> <span class="n">latmax</span><span class="p">)</span>
<span class="n">dx</span> <span class="o">=</span> <span class="p">(</span><span class="n">xmax</span><span class="o">-</span><span class="n">xmin</span><span class="p">)</span> <span class="o">/</span> <span class="n">maxn</span>
<span class="n">dy</span> <span class="o">=</span> <span class="p">(</span><span class="n">ymax</span><span class="o">-</span><span class="n">ymin</span><span class="p">)</span> <span class="o">/</span> <span class="n">maxn</span>
<span class="n">dx</span> <span class="o">=</span> <span class="nb">max</span><span class="p">(</span><span class="n">roaring</span><span class="o">.</span><span class="n">dx</span><span class="p">,</span> <span class="n">dx</span><span class="p">)</span>
<span class="n">dy</span> <span class="o">=</span> <span class="nb">max</span><span class="p">(</span><span class="n">roaring</span><span class="o">.</span><span class="n">dy</span><span class="p">,</span> <span class="n">dy</span><span class="p">)</span>
<span class="n">x</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">arange</span><span class="p">(</span><span class="n">xmin</span><span class="p">,</span> <span class="n">xmax</span><span class="p">,</span> <span class="n">dx</span><span class="p">)</span>
<span class="n">y</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">arange</span><span class="p">(</span><span class="n">ymin</span><span class="p">,</span> <span class="n">ymax</span><span class="p">,</span> <span class="n">dy</span><span class="p">)</span>
<span class="n">yy</span><span class="p">,</span> <span class="n">xx</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">meshgrid</span><span class="p">(</span><span class="n">y</span><span class="p">,</span> <span class="n">x</span><span class="p">)</span>
<span class="n">lons</span><span class="p">,</span> <span class="n">lats</span> <span class="o">=</span> <span class="n">transformer</span><span class="o">.</span><span class="n">transform</span><span class="p">(</span><span class="n">xx</span><span class="p">,</span> <span class="n">yy</span><span class="p">,</span> <span class="n">direction</span><span class="o">=</span><span class="s1">'inverse'</span><span class="p">)</span>
<span class="n">extent</span><span class="o">=</span><span class="p">[</span><span class="n">xmin</span><span class="p">,</span> <span class="n">xmax</span><span class="p">,</span> <span class="n">ymin</span><span class="p">,</span> <span class="n">ymax</span><span class="p">]</span>
<span class="n">transform</span> <span class="o">=</span> <span class="n">crs_plot</span>
<span class="k">else</span><span class="p">:</span>
<span class="n">dlon</span> <span class="o">=</span> <span class="p">(</span><span class="n">lonmax</span> <span class="o">-</span> <span class="n">lonmin</span><span class="p">)</span> <span class="o">/</span> <span class="n">maxn</span>
<span class="n">dlat</span> <span class="o">=</span> <span class="p">(</span><span class="n">latmax</span> <span class="o">-</span> <span class="n">latmin</span><span class="p">)</span> <span class="o">/</span> <span class="n">maxn</span>
<span class="n">dlon</span> <span class="o">=</span> <span class="nb">max</span><span class="p">(</span><span class="n">roaring</span><span class="o">.</span><span class="n">dx</span><span class="p">,</span> <span class="n">dlon</span><span class="p">)</span>
<span class="n">dlat</span> <span class="o">=</span> <span class="nb">max</span><span class="p">(</span><span class="n">roaring</span><span class="o">.</span><span class="n">dy</span><span class="p">,</span> <span class="n">dlat</span><span class="p">)</span>
<span class="n">lons</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">arange</span><span class="p">(</span><span class="n">lonmin</span><span class="p">,</span> <span class="n">lonmax</span><span class="p">,</span> <span class="n">dlon</span><span class="p">)</span>
<span class="n">lats</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">arange</span><span class="p">(</span><span class="n">latmin</span><span class="p">,</span> <span class="n">latmax</span><span class="p">,</span> <span class="n">dlat</span><span class="p">)</span>
<span class="n">lats</span><span class="p">,</span> <span class="n">lons</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">meshgrid</span><span class="p">(</span><span class="n">lats</span><span class="p">,</span> <span class="n">lons</span><span class="p">)</span>
<span class="n">extent</span><span class="o">=</span><span class="p">[</span><span class="n">lonmin</span><span class="p">,</span> <span class="n">lonmax</span><span class="p">,</span> <span class="n">latmin</span><span class="p">,</span> <span class="n">latmax</span><span class="p">]</span>
<span class="n">transform</span> <span class="o">=</span> <span class="n">crs_lonlat</span>
<span class="n">img</span> <span class="o">=</span> <span class="n">roaring</span><span class="o">.</span><span class="n">mask</span><span class="o">.</span><span class="n">contains_many</span><span class="p">(</span><span class="n">lons</span><span class="o">.</span><span class="n">ravel</span><span class="p">(),</span> <span class="n">lats</span><span class="o">.</span><span class="n">ravel</span><span class="p">())</span><span class="o">.</span><span class="n">reshape</span><span class="p">(</span><span class="n">lons</span><span class="o">.</span><span class="n">shape</span><span class="p">)</span><span class="o">.</span><span class="n">T</span>
<span class="kn">from</span><span class="w"> </span><span class="nn">matplotlib</span><span class="w"> </span><span class="kn">import</span> <span class="n">colors</span>
<span class="n">cmap</span> <span class="o">=</span> <span class="n">colors</span><span class="o">.</span><span class="n">ListedColormap</span><span class="p">([</span><span class="n">ocean_color</span><span class="p">,</span> <span class="n">land_color</span><span class="p">])</span>
<span class="n">ax</span><span class="o">.</span><span class="n">imshow</span><span class="p">(</span><span class="n">img</span><span class="p">,</span> <span class="n">origin</span> <span class="o">=</span> <span class="s1">'lower'</span><span class="p">,</span>
<span class="n">extent</span><span class="o">=</span><span class="n">extent</span><span class="p">,</span> <span class="n">zorder</span><span class="o">=</span><span class="mi">0</span><span class="p">,</span> <span class="n">cmap</span><span class="o">=</span><span class="n">cmap</span><span class="p">,</span>
<span class="n">transform</span><span class="o">=</span><span class="n">transform</span><span class="p">)</span>
<span class="k">if</span> <span class="n">fast</span><span class="p">:</span>
<span class="n">show_landmask_roaring</span><span class="p">(</span><span class="n">get_mask</span><span class="p">())</span>
<span class="k">else</span><span class="p">:</span>
<span class="k">if</span> <span class="n">latmax</span><span class="o">-</span><span class="n">latmin</span> <span class="o"><</span> <span class="mi">3</span> <span class="ow">and</span> <span class="n">lonmax</span><span class="o">-</span><span class="n">lonmin</span> <span class="o"><</span> <span class="mi">3</span><span class="p">:</span>
<span class="c1"># Check first if there is land at all within map bounds</span>
<span class="n">roaring</span> <span class="o">=</span> <span class="n">get_mask</span><span class="p">()</span>
<span class="n">maxn</span> <span class="o">=</span> <span class="mf">512.</span>
<span class="n">dlon</span> <span class="o">=</span> <span class="p">(</span><span class="n">lonmax</span> <span class="o">-</span> <span class="n">lonmin</span><span class="p">)</span> <span class="o">/</span> <span class="n">maxn</span>
<span class="n">dlat</span> <span class="o">=</span> <span class="p">(</span><span class="n">latmax</span> <span class="o">-</span> <span class="n">latmin</span><span class="p">)</span> <span class="o">/</span> <span class="n">maxn</span>
<span class="n">dlon</span> <span class="o">=</span> <span class="nb">max</span><span class="p">(</span><span class="n">roaring</span><span class="o">.</span><span class="n">dx</span><span class="p">,</span> <span class="n">dlon</span><span class="p">)</span>
<span class="n">dlat</span> <span class="o">=</span> <span class="nb">max</span><span class="p">(</span><span class="n">roaring</span><span class="o">.</span><span class="n">dy</span><span class="p">,</span> <span class="n">dlat</span><span class="p">)</span>
<span class="n">lons</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">arange</span><span class="p">(</span><span class="n">lonmin</span><span class="p">,</span> <span class="n">lonmax</span><span class="p">,</span> <span class="n">dlon</span><span class="p">)</span>
<span class="n">lats</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">arange</span><span class="p">(</span><span class="n">latmin</span><span class="p">,</span> <span class="n">latmax</span><span class="p">,</span> <span class="n">dlat</span><span class="p">)</span>
<span class="n">lats</span><span class="p">,</span> <span class="n">lons</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">meshgrid</span><span class="p">(</span><span class="n">lats</span><span class="p">,</span> <span class="n">lons</span><span class="p">)</span>
<span class="n">land</span> <span class="o">=</span> <span class="n">roaring</span><span class="o">.</span><span class="n">mask</span><span class="o">.</span><span class="n">contains_many</span><span class="p">(</span><span class="n">lons</span><span class="o">.</span><span class="n">ravel</span><span class="p">(),</span> <span class="n">lats</span><span class="o">.</span><span class="n">ravel</span><span class="p">())</span>
<span class="k">if</span> <span class="ow">not</span> <span class="n">land</span><span class="o">.</span><span class="n">any</span><span class="p">():</span>
<span class="n">logger</span><span class="o">.</span><span class="n">debug</span><span class="p">(</span><span class="s1">'No raster land within plot area, skipping plotting of GSHHG vectors'</span><span class="p">)</span>
<span class="k">return</span>
<span class="n">land</span> <span class="o">=</span> <span class="n">LandmaskFeature</span><span class="p">(</span><span class="n">scale</span><span class="o">=</span><span class="n">lscale</span><span class="p">,</span> <span class="n">facecolor</span><span class="o">=</span><span class="n">land_color</span><span class="p">,</span> <span class="n">globe</span><span class="o">=</span><span class="n">crs_lonlat</span><span class="o">.</span><span class="n">globe</span><span class="p">,</span> <span class="n">levels</span><span class="o">=</span><span class="p">[</span><span class="mi">1</span><span class="p">,</span><span class="mi">5</span><span class="p">,</span><span class="mi">6</span><span class="p">])</span>
<span class="n">ax</span><span class="o">.</span><span class="n">add_feature</span><span class="p">(</span><span class="n">land</span><span class="p">,</span> <span class="n">zorder</span><span class="o">=</span><span class="mi">0</span><span class="p">,</span>
<span class="n">facecolor</span><span class="o">=</span><span class="n">land_color</span><span class="p">,</span>
<span class="n">edgecolor</span><span class="o">=</span><span class="s1">'black'</span><span class="p">)</span></div>
<div class="viewcode-block" id="Reader">
<a class="viewcode-back" href="../../../autoapi/opendrift/readers/reader_global_landmask/index.html#opendrift.readers.reader_global_landmask.Reader">[docs]</a>
<span class="k">class</span><span class="w"> </span><span class="nc">Reader</span><span class="p">(</span><span class="n">BaseReader</span><span class="p">,</span> <span class="n">ContinuousReader</span><span class="p">):</span>
<span class="w"> </span><span class="sd">"""</span>
<span class="sd"> The global landmask reader is based on the coastline data from</span>
<span class="sd"> GSHHG (https://www.ngdc.noaa.gov/mgg/shorelines/) optimized for</span>
<span class="sd"> checking against landmasks.</span>
<span class="sd"> """</span>
<span class="n">name</span> <span class="o">=</span> <span class="s1">'global_landmask'</span>
<span class="n">variables</span> <span class="o">=</span> <span class="p">[</span><span class="s1">'land_binary_mask'</span><span class="p">]</span>
<span class="n">proj4</span> <span class="o">=</span> <span class="kc">None</span>
<span class="n">crs</span> <span class="o">=</span> <span class="kc">None</span>
<span class="k">def</span><span class="w"> </span><span class="fm">__init__</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
<span class="bp">self</span><span class="o">.</span><span class="n">proj4</span> <span class="o">=</span> <span class="s1">'+proj=lonlat +ellps=WGS84'</span>
<span class="bp">self</span><span class="o">.</span><span class="n">crs</span> <span class="o">=</span> <span class="n">pyproj</span><span class="o">.</span><span class="n">CRS</span><span class="p">(</span><span class="bp">self</span><span class="o">.</span><span class="n">proj4</span><span class="p">)</span>
<span class="nb">super</span><span class="p">(</span><span class="n">Reader</span><span class="p">,</span> <span class="bp">self</span><span class="p">)</span><span class="o">.</span><span class="fm">__init__</span><span class="p">()</span>
<span class="c1"># Depth</span>
<span class="bp">self</span><span class="o">.</span><span class="n">z</span> <span class="o">=</span> <span class="kc">None</span>
<span class="bp">self</span><span class="o">.</span><span class="n">xmin</span><span class="p">,</span> <span class="bp">self</span><span class="o">.</span><span class="n">ymin</span> <span class="o">=</span> <span class="o">-</span><span class="mi">180</span><span class="p">,</span> <span class="o">-</span><span class="mi">90</span>
<span class="bp">self</span><span class="o">.</span><span class="n">xmax</span><span class="p">,</span> <span class="bp">self</span><span class="o">.</span><span class="n">ymax</span> <span class="o">=</span> <span class="mi">180</span><span class="p">,</span> <span class="mi">90</span>
<span class="c1"># setup landmask</span>
<span class="bp">self</span><span class="o">.</span><span class="n">mask</span> <span class="o">=</span> <span class="n">get_mask</span><span class="p">()</span>
<div class="viewcode-block" id="Reader.__on_land__">
<a class="viewcode-back" href="../../../autoapi/opendrift/readers/reader_global_landmask/index.html#opendrift.readers.reader_global_landmask.Reader.__on_land__">[docs]</a>
<span class="k">def</span><span class="w"> </span><span class="nf">__on_land__</span><span class="p">(</span><span class="bp">self</span><span class="p">,</span> <span class="n">x</span><span class="p">,</span> <span class="n">y</span><span class="p">):</span>
<span class="n">x</span> <span class="o">=</span> <span class="bp">self</span><span class="o">.</span><span class="n">modulate_longitude</span><span class="p">(</span><span class="n">x</span><span class="p">)</span>
<span class="n">x</span> <span class="o">=</span> <span class="n">x</span><span class="o">.</span><span class="n">astype</span><span class="p">(</span><span class="n">np</span><span class="o">.</span><span class="n">float64</span><span class="p">)</span>
<span class="n">y</span> <span class="o">=</span> <span class="n">y</span><span class="o">.</span><span class="n">astype</span><span class="p">(</span><span class="n">np</span><span class="o">.</span><span class="n">float64</span><span class="p">)</span>
<span class="k">return</span> <span class="bp">self</span><span class="o">.</span><span class="n">mask</span><span class="o">.</span><span class="n">contains_many</span><span class="p">(</span><span class="n">x</span><span class="p">,</span> <span class="n">y</span><span class="p">)</span></div>
<div class="viewcode-block" id="Reader.get_variables">
<a class="viewcode-back" href="../../../autoapi/opendrift/readers/reader_global_landmask/index.html#opendrift.readers.reader_global_landmask.Reader.get_variables">[docs]</a>
<span class="k">def</span><span class="w"> </span><span class="nf">get_variables</span><span class="p">(</span><span class="bp">self</span><span class="p">,</span>
<span class="n">requestedVariables</span><span class="p">,</span>
<span class="n">time</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span>
<span class="n">x</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span>
<span class="n">y</span><span class="o">=</span><span class="kc">None</span><span class="p">,</span>
<span class="n">z</span><span class="o">=</span><span class="kc">None</span><span class="p">):</span>
<span class="w"> </span><span class="sd">"""</span>
<span class="sd"> Get binary mask of whether elements are on land or not.</span>
<span class="sd"> Args:</span>
<span class="sd"> x (deg[]): longitude (decimal degrees)</span>
<span class="sd"> y (deg[]): latitude (decimal degrees)</span>
<span class="sd"> ...</span>
<span class="sd"> x, y is given in reader local projection.</span>
<span class="sd"> Returns:</span>
<span class="sd"> Binary mask of point x, y on land.</span>
<span class="sd"> """</span>
<span class="bp">self</span><span class="o">.</span><span class="n">check_arguments</span><span class="p">(</span><span class="n">requestedVariables</span><span class="p">,</span> <span class="n">time</span><span class="p">,</span> <span class="n">x</span><span class="p">,</span> <span class="n">y</span><span class="p">,</span> <span class="n">z</span><span class="p">)</span>
<span class="k">return</span> <span class="p">{</span><span class="s1">'land_binary_mask'</span><span class="p">:</span> <span class="bp">self</span><span class="o">.</span><span class="n">__on_land__</span><span class="p">(</span><span class="n">x</span><span class="p">,</span> <span class="n">y</span><span class="p">)}</span></div>
</div>
</pre></div>
</div>
</div>
<footer>
<hr/>
<div role="contentinfo">
<p>© Copyright 2020, Knut-Frode Dagestad (knutfd@met.no) and Gaute Hope (gauteh@met.no)..</p>
</div>
Built with <a href="https://www.sphinx-doc.org/">Sphinx</a> using a
<a href="https://github.com/readthedocs/sphinx_rtd_theme">theme</a>
provided by <a href="https://readthedocs.org">Read the Docs</a>.
</footer>
</div>
</div>
</section>
</div>
<script>
jQuery(function () {
SphinxRtdTheme.Navigation.enable(true);
});
</script>
</body>
</html>