{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Intersecting SIDs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import matplotlib.tri as tri\n",
    "import cartopy.crs as ccrs\n",
    "import numpy\n",
    "import pystare"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false,
    "jupyter": {
     "outputs_hidden": false
    },
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "proj = ccrs.PlateCarree()\n",
    "transf = ccrs.PlateCarree()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false,
    "jupyter": {
     "outputs_hidden": false
    },
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "def plot1(lon, lat, lons, lats, triang, c0='r', c1='b', transf=None, lw=1, ax=None):\n",
    "    if lon is not None and ax is not None:\n",
    "        x = numpy.zeros([lon.size + 1], dtype=numpy.double)\n",
    "        x[:-1] = lon[:]\n",
    "        x[-1] = lon[0]\n",
    "        y = numpy.zeros([lat.size + 1], dtype=numpy.double)\n",
    "        y[:-1] = lat[:]\n",
    "        y[-1] = lat[0]\n",
    "        ax.plot(x, y, True, transform=transf, c=c0)\n",
    "    ax.triplot(triang, c1 + '-', transform=transf, lw=lw, markersize=3)\n",
    "    ax.scatter(lons, lats, s=10, c=c1, transform=ccrs.PlateCarree())\n",
    "    return"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "resolution = 6\n",
    "resolution0 = resolution\n",
    "lat0 = numpy.array([10, 5, 60, 70], dtype=numpy.double)\n",
    "lon0 = numpy.array([-30, -20, 60, 10], dtype=numpy.double)\n",
    "hull0 = pystare.cover_from_hull(lat0, lon0, resolution0)\n",
    "\n",
    "lons0, lats0, intmat0 = pystare.triangulate_indices(hull0)\n",
    "triang0 = tri.Triangulation(lons0, lats0, intmat0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "resolution1 = 6\n",
    "lat1 = numpy.array([10, 20, 30, 20], dtype=numpy.double)\n",
    "lon1 = numpy.array([-60, 60, 60, -60], dtype=numpy.double)\n",
    "hull1 = pystare.cover_from_hull(lat1, lon1, resolution1)\n",
    "\n",
    "lons1, lats1, intmat1 = pystare.triangulate_indices(hull1)\n",
    "triang1 = tri.Triangulation(lons1, lats1, intmat1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false,
    "jupyter": {
     "outputs_hidden": false
    },
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "r01 n01:  79\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAGgAAADnCAYAAADy+hOsAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAA0HklEQVR4nO1dd3hUxfp+d9N7JyRBSoAUihI6MUBAKZFQJAmCoj8UrCggXCwRFQXxAiIXVBQNV3oLLZQLioQrIj3ghdACQgBpCSEhlPR9f3/Mns227G6yu8mqeZ9nnt1zzpw5c+Y9880338x8IyOJetgu5HWdgXoYRj1BNo56gmwc9QTZOOoJsnHYmxKpf//+vHXrlsa5iooK5Obm4ubNmwgNDYWHh4dFMlRUVIRTp06pjiMiIuDm5maRtGsLJSUlyMzMBAA4ODhAJpMhICAAXl5ecHFxUcXLyMj4gWR/g4mRNBo6dOjAqhAUFEQA/Oabb1hUVFRlPFNw4cIFAlCF3Nxcs9IjyYYNSYCUy8WvTEaGhJDe3mRiYmW8xETS2zufvXod4C+//MKYmB309s5lYqK45uurGb8qnDt3jgDo4eHBAQMGsEGDBrx165ZOvPz8fAI4QiNlbxJBzZs3rzJDqampGoV67do1E4pNP9TTEd9OzfD222SbNuLtDIdCAl8TGEmgPwF/Ah0IdCPQSy0vGwjcIlBhlLCioiK+9dZbbN68uer+qVOnMiMjg7Nnz+a6det48uRJzpw503IEAWB6ejofPHigN0P9+vVTZebIkSMmF+TRo0e5ePFizpw5ky+99JIqDZlMxszMzGpQQiYnC1Kio00h5iqBtwj4EkgikEJgG4EsrXifE4hVBh8CjgR6EphN4DKB32lv34rOzsns2/ck4+Nv0N09lS1ajCMABgUF8dVXX2WvXpVkDxkyhM2aNZOOLUOQs7Mzo6KiGBoayjFjxvCPP/7QKJwpU6YQAPfv329SYWZlZfGVV15hSEgIR44cyUmTJnH+/Pn83//+x/Lych45coSZmZl88OABL126xCNHjjA1NVWvCE1OJv39TSFFCtsJBBB4g8AFrWuKKu75lsAHBLYQWERA+iDf1ar1TgQGEPgngbYMDOzH+fPnc/Pmzfziiy94+vRpVb5PnTpl2Rrk6urKhx9+mAB44cIFjUKaP38+AWhkQB0KhYKzZ89mdHQ0vby8GBQUxHHjxjE/P18nrnpNcnNzY0BAAFu2bEkAbNmyJcvKyjTIMZ2YCgLTCAQR2KN17X8E2iuf25ZAqdb1wVpESKEHATclIVMJXNR5roODaPP0iUSLEeTg4EAAjI2NZXZ2tl4SJk+ezMTERI0ClPCf//yHPj4+TE9PZ25uLhUKhd40SM027ZVXXuHMmTOZlJRENzc3hoaGatQiR0dTiHlAYD1F29KDwFWV4gCUE9hIwFX5TDmBYQTK9KTTR4ucydX4OCpDRIQVCJIytXjx4ioL9s6dO+zbty9jYmI0CnHx4sUEwLS0tCrv1UZWVhZnz57Njz76iBMmTOCSJUs0am1yMunnZ0qBfErAi6INWausReJajx4KymR9CXQkEEqgM4HrBETaJJmQoJ5WKYF9BC7ViBj1kJxsYYL8/f2ZkpLCEydOGCzYiooKjhgxgpGRkezXrx+//vprPv/88xw7dqzJ5BhD9cRaMIEM1bGTE3n9ulTzcgjYEdhL4FPKZGEEjlMuX0fgn/Tw+BeBe3RyqlTRLRXatLEwQS1btqSbmxudnZ157NgxgwWoUCj4yy+/cOPGjYyPj2fv3r11lApzYG+v/cKFBP5D4JCewuhJIcKuEBhNYLratQIlgT0plIWe1G1jfrQoMVKQ2iKLanGpqalct24dAwMDq60CWwqJifpe+HG1Al2jdW0thSpdWehbt/7Ks2fJ0NDpSlKuEsiiXD5AGectAqcJ3LQKOQApCRSLEaRuSUhISKCbmxt///332mFFieRkfaLmjrJQn1X+jtYRJX5+ZRokBQYGcvz48QwMDKS/fzoBUi4fRuApCvX5tNWIAUhXV1Jqjq1C0K1btzhlyhT6+PiwadOm1eqY1hQREdov+kD5pbur1Y4QJRGpBEhPT3HvnTt3CDhTW3zZ27cmEKEMoND0LEeEXC5Ua0nFjo4WNUddV7IKQRIUCgXnz59PJycnZmVlWYcZViXW9hLwpOjPSIXuRuAn5fnP6OV1mY888oja9U84fPhwTp48mU5OfSmsAQcJLCXwR41IsLMjXVzEB2SqrU4dViNIoVBw7dq1fOyxxwiAy5cvrzkDRl9CX1ikLHRvSp1IO7tnKNTgAAKgt3d7vvPOOywvL2d2trjv/HnxRZtbO9T7Mua9mwWVBOlLdHFxYYMGDdiqVSvVuT59+lRpRTAH+msPKfozVwgUUSbrqsqHXN6AMlkAnZwacMWKFRqdZkuJLakPYwlYvKPq5uZGANy8eTMBMDQ0lF26dFEVkD4rgjlwcDBeaPb22VRvW7755hvm5t7ijRvkvn3k8uX6VHPTg7u7rvgyZAmpDixGUPv27ZmQkEAXFxcOHTqUTZs2ZZMmTdixY0eNwrFkf6fq2qMdFBR2sMYEJpgQv4zqRtEe2K1x3de30jL+9NPpfOqpp7hmzRpu376dCQkJlOyS2lAoFFy1ahX79+/Phx9+mJMnTzb6jhYjyMnJSYMIAJw7dy5J8uLFizx8+DB37txpJiWakMlM/8plqOAs/B8dIOMfcNOJUABwLp7j72r5b4FUXkUQFQCLYM9LCOEtWaXatbhjR3p6erJRo0Y6775161aNvBYXF7NNmzZs0aIFO3XqRABs1KiR0Xe0uIgDwI8++ogZGRkWokE/EhPJBkq7mHaIwxaWAfwUYC+AQwGOAfgCQDeAFcqQAvAdZfAD6ADwM62C3qCMq/2Q28q0vgNYbG/PeGdnxtrbc3TTpgwMDKSdnR2XLFnCiooKkoIgAMzKylKlffbsWaPvaTGC2rVrxy1btrCkpMTCVKghMVE0Og4OLAdYCnt2xgHKUMF1GKpRgPeUhbAUYCrArwEuAHgAYDnAJQDdAX4M8E2ATynjN1AeuyiPgwGGAbymRdAygD4AMwH+B+AogIOUH8UlgF/J5fRycGBKSooq+2PGjOHgwYN58uRJFUlnzpwx+MpW7QdZDImJ1ZNnyhCuJET9XBnAvsoak6p2XiqwnwB206pFALgS4Go95wGwHcCvlOmNBugP0ANgBMCDgPioEhNZVFTExo0b09vbm/Hx8QTArl27Gnx12yfI2bnaxBBglrLw7mmdTwPoCvAfmE4CHINvSYCdISMAFgI8A3C+GgG9AD6kLHCo1bYNAIuqma9eAGPlcm6KjiYAo4ZlUwiq/Xlx770HtG0L2NsDxcU1SqKh8jcBwC8ADgFYhAYYDyAS7RCPKRiFfkiBHMPRGDEgugIYByBC+QsA/wWQDiAWwBnluaHK33IAtwD8CkAB4GcAp43kqxCATKFA2b59cAFwa9w4I3eYAGMM0lI1qIaiTApn0JIP4KRsI57mAoBP6xFJvgC7AIxGMJ8C+D7A8QC9ldd/AfiSXM6Vrq5kbCwJ8BxAOcBGAHsDXAiwszJ+QwgFQ0rfUB5vA5wDcAjAlwHeBMQ7V2EDQp2LuMTEGouxUJznXmhO0bkGMVZdOiSJZxHEt5QFPxbgCIDTAZ6tIr0PAQYZeF6JkiQAnKk8VwxwsfJcE4CtavhxERDqu00RZHpPU2/ohl9JgIvwPO8ra44UHkAoCQkQWtVzWKw3DYXa/58BekJobQkA9yvP50HUxBglMdla90lEPaYky5x30q5NdUOQ7thAtcJofMcbaGAwzodow55qxx9jiklplwI8BnAehBo9CkIcdQe4WUmcHEJFj4NQHgAwWvn7wEj6fwC8b8p7Kq2ttU+Qu3uNiZG+2m2IMxpvMJozEWAAbvI6Amv0vHNKonoBTNe6dgvgKgg1+gOAXQG+BlHbqkpvGdSmLJuSh+TkWiRI/3CnxcMegFEA2wO8aOVnSSFdreDXGoh3GGAghLh8DOBJiH7Z5aru8fSsJYLMFGmmhmMQloAlEG1CbTxTCq2UBAVB1K6q4r2ujPcYhAiViP0aYL5W3AdyuZUJSk4mPTxqpYAUAPtBaGna1wZgC3Ng0iS5GoffIawTsQCHQ7/9jhAdW4mU79T+R0LY9hoqiXMGaAcLzs3WIaiWao0UeuIZ+iOEQD79kEs5yjWiaGtdBITIlebdKm18Gv8jIsSYgjRerX5NT9dgDIR4bQth46sqrz9BmIJ6Q5iQXoKwbswBeAWiBt4DmDthgpUIMlN9NjWcRzPOwiQCCrrDi8AmvVHz4S6GCuBbWeiWGPaU1pgoycoG2BLgKwBbV/NdflTWJJXpSCYjSSsRZIY1oKowCv+mE4r0Xm6M8/SGC8XKgqqTqe6EjZqQ9aGHB5+XydgE4HQTy+EyhEibp35e2Wm1GEHu7u6cNm1adefdmhSK4GRCtGkEFlR53VKTOIxh69atjIiI4NGjRxkYGMijR48aNGFdhRCJMdrX3N1NJsgkYylJ3LlzB/jii2rb+oxhKj40IVaeMugiIgI4bcyKaSHExcUhJiYG8+bNw4wZMzBkyBAULV0KKJStYHIy0KYNzrm6YjKAlgCaAfhJO6HAQNMfakoNUok4OzuL1p5S2DPM6EzO+wRaEPhB55pSlNcqcnJy6O3tzfXr19PZ2Zm7d+9WXSsoKOCiRYvo5+fHV0JDudPBQbd/KJerZi/C4m2QpchxcSGTk01c37OPACgmyWtes2q7YwAHDhxQLWP88ssvuW3bNk6cOJENGzbkkCFDmJqaqjvzJy1NZ2qp7RBkb6+jWVWvOWtOYKdOknWJK1euKD8csEWLFnzrrbe4YcMGpqSkmDwd2mIE+fj4cOrUqRxsZ0dAWJJrSowE01ZhS+E7ikW8m1XnLDmBsKYYPnw4ATAxMZFNmzaln58fAfD999836X6LEdSoUSMOGzaML3h7EwBlEJ2umhAjofoK4Y8Uy+RPSEpQnaOsrIzz5s3jkiVLeOLECSoUCnbu3Jk7duww6X7Lizhlqb4EMXK5Q7sUAwKq9WlX376aTCCRkZGFJj+jNpGRkcGQkBCWlpaaFN/yBJEqRwTfQ5gwhgMcI5NxQ0JCtafEVt8oUUjgKfr6hqkmZJSUlDA7O5slJSW8desWP/vsM/bu3ZutW7fmmDFjVHPXrA2FQsH+/ftz3rx5Jt9jHYJIUUsCAnjE25vz+/Thv/71Lz7yyCNMSEhgeXl5tV6suma95GRy6dKl9PPzY3x8PL28vBgcHEwHBwd6enryueee49atW3n06FFGRUUxPj6eTZs25dq1a6uVr+pi4cKFbNOmTbXmDlqPID24e/cuw8LCuG3btmq8lkByslhwZYwcdbX68uXLXLp0qcr1zP3795mXl6eRbnZ2NufMmcMPP/yQLi4uLCy0jmjctGkTGzRoYHSiojZqlSCS/Pzzz9msWbNqZVIbycn6VyOYY87Jycmho6OjxVdfFBUVcezYsWzUqBEPHjxY7ftNIcii8+LkcjkuXryIxx9/HDt37kRZWVm10/jkE6CsTFhNAgJESE42z5xz/PhxhIWFwd7eJO9rJqG0tBQDBw7EjRs3cPz4cXTu3NliaatDJog0EkkmO2KVp9fjFo34izOJoHrUHepdYto46gmycdQTZOOoJ8jGUU+QjaOeIBtHPUE2DoNda30Ozf8qIImjR48CADp06GDyfRcuAIWFgKcnEBqqe6wPhYWFOHfuHORyORQKBYKDg+Hk5ARfX18ARhycG7IDWXXxsA0AyiHr6OgtGt5EtP1ix8eTPj5kUJCmfTAgQPO4Tx/ywgUyJ4ccMkQYgBMTyffee48AmJeXRz+/EMpkYmR68ODBUj6qtMn9LQjSV+De3goCYHj4UgJ3VYXspLlWrBphFWWyf7BBgx/p4KDp3tnR8SM6O7dg584kcJTACNXHsWjRor8+QRUVFTx06BCPHz9OUhDh41PMQYPInj01C7Jy5tg5Ao1qSIYUSimmg71IIJJy+VQKt87NKFyjiXhubv+jl5cfvb3fIZCvPL9RRdJfjiCpRgwZcp8pKSl0cwsn4E4Hh0DKZG4UHuLlFD5J+yi/2C4UTseTCHxE4UpztJ5CVw8KCkd/EQQ6Ke+fSeHv9DCBRyi8Bb9DIFftni8o3G2KmiSTkU5OFyiXv0DAj0AWXVwUHDlyhuS1/89NkLqI0nSV/JbyC5xJoEj5RedTeGQspXA0vo3Cad9uAssI/JvA2xRTie8YIWgJxaTJDcov/kUKt5kyAg0IfEVNb/U/Ubh53q7M19c6aT700DRGRr6gmh6nUCj+fASpEzJokKEC3KEsiOU1FFHGwlQCMXrOV1B3G4HxrBRZ0+nomKD834ZAniqet/dNNmzYkPPmzVPNl6gzgvRMpjRKSHx8dQsRFHsrWJqcCmVNAyUxpX5dd758BIX48qa0RYC7+yHa2Y1Q1sLfCIh3PHfuHKOjoxkeHs7+/fvXDUFpaZUO+VxcNElSJ2Tw4OoXnqajP1ByJGt+UBA4S2AxgVYE2lG9sTd831wCoL39wzrXg4NXUibzY7NmL3LhwoW8fv06FQoFDxw4wNdee61uCBo7VjOTcrmYsaXZhpgWDM/Zn0Qhiswlp5RCsfAnMJzC53ZVu6GIJT7h4dLxFAIPUyZzIQA2bKj5wbRpIyawvP322xw0aBAdHBw4c+ZMFhQUsLy8vG4I6tLFtILRN3mxehMau1ugBikIPE9RG6/puV5EoFxDrMnlQsx5epJ2dssICB9xEydOZNu2/6dxv/pczqysLJWnyh49eki+5qokyCq2uKQk4OBB0+IqFKad0498AIcBRJt6QxV4B8BJCHdIQcpzxwD4A/AE4AIgHKRmHklh5hk0KAkA4OzsjOLiYiQlNZeWCiE5WUyEkdCyZUscPnwY165dQ2lpKby8vAxnzRo1yNfXXHFjaniOwFCDosh4uEKxPY2kaeUSuE0gTFmjQOAXAj8xKIi8cUO3hvv4PCAArlq1im3atOGyZctMLquFCxeaL+KuXr3K1NRUbtu2jR9//DFPnjxp8KFPPFEb5BSpFaD2TlrVCWUUnuk/ojD5QEl8KYUav18jvp2drveB+PgcSq6ro6Ojq7WP340bN8wnSH3DvIEDBxqcPZqYKLQ2V1drE1SuRtB8M9JZoUxjmvJ4AqtaUa4etHfVysvL09mZTB0nTpzgihUruH//fubl5fHw4cOcO3cuw8LCzCdI2ldOfZPBgoIChoSEMDIyUrXVWS2t0FcLxRTrhv5bg3sVBI4rydHeMk03aBtRfX1NqiAkVVtyEhDegD08PNiiRQsCoJ1Yc2WekmBnZ4cOHTpobNLq6emJZ599Fp07d4arqysAID3dlNQsCUcADyB8Lp4AsADAWwDuGbmvAsAQAI8BaApTlAztSbKdO/+Oe/fu4eTJk6rC3L59O+Li4vDQQw8hOzsbFy9eBAC4uroiKSkJ3t7eiI6Oxrp165CZmQkHBwesW7fO8INNqUGGkJyczKVLl9ZRDZK2p+lGsfvJcxR9GRDIVNaSCxT9pHEU5pgXKexl0dTdq05Bsf1mBoU97S0K42oEgRZ0cOhIZ+f36e/fU1UjoKwV3t7eDA8P5+uvv65xbf/+/UxLS2NJSQkLCwv56aefMioqih06dCAALliwwHwRZwhlZWWqdUFSJ7QWHF8ZCFLb9C6FBdtbScrnFEbV2RRG0xta93yhJCKQwkr9OO3tJ1K0UUcJZNHd/T98/fXXuXnzZt69e5f5+fksLS3l8ePHmZOToyqHvLw8AqK9Lisr01lUnJ+fz/nz5/PDDz9kWVmZdQmSUPu1x1A4qPYV5xiId4rACxS2so4UbVGlyq5tb6vOqvImTZpw+vTpeq/985//ZIsWLVSk1QpBtdf3MSXMIdCbYvvnKALvUezS1Y9iXGcZhVrtR+ATiv2I7lWZntI1drUg7Uuub6eykpISjQ0bDRFkMUtCr16WSslclACYD2AGgP8AmALADkA3AH0AhCvPXwNwAEAygEcBuKlS0F6l4uEBpKZWLxd79+6Fn58fDhw4oHPN0dERbdq0MSkdiy2YKSiwVErmYgaAdgC6KI+HotIbds3Qu3f14k+aNAkbN27EZ599hnbt2pn1bIuIuLpvf8oplukPIdCEwnxjXprBwTXb9oykaosaPz8/3r9/32h8WLsNqtv25zCBphSN/FcUvn3MT7c6HVFt5ObmSnMNGBISwp07d1KhULCsrEzvqnOrExRYM8e7FFqU/rF708KvFHMD1plNiLalwBJ+gC5fvsxp06YxMjKSwcHBdHV1ZUBAAPv27cvU1FTm5uZanyDzxFspK7Womtw/gMCXFqkxcrn+7dAsgYqKCl68eJEFBQW8dOkSV69ezZ49e9Ld3Z0rVqywLkF1J95SKeaf5Vs03dr0oLV06VI+/vjjBgkyW83u3t3cFGqKhQA+B+Bt0VRr056YlJSEiooKg3HMJignB5DJzE2lJsiF6N+YBzutJKqrUpsDZ2dn7Nq1y2AcswhKSgL2K3epqH30A7Dd7FSkD9jZGUhMrH6H1FzIjHzdZhFkhHwrIgfAGgADapyCk5PmcYsWtU+OKTCLoGosq9GDcgCGPJHshpiskaN1PgfCMvAMqkuQj48gJjoamDRJ89qgQdVKqtZglqnn5k3R/tRMxL0KMdC2D8AmiIJ/U3mtFIKAIGgqARsBvAhgDIBp1X5ifr6uGNu8WZCjPvPGplCVekcjarb55p3tyj7MIwSg/O1CMTzwhPLczxROZWdQWKeDKQbTav5ccywE1gKsoWbv3GnOZ1EBIA5Cwr4C4FOI2nMQwFkASQB2QgxpD4TY7u81ABcAtDfnwbWqpVkCNRJxSUnA3bvmPNYOYt/FUAAOaudHQJASozyOBzAdgkTz4OwMxMfbpiJgCNUmKCkJMDbPwTSE6znnqPy9C2AugCMQ2lr1YWdXqUIDtqulGUO1CHrvPWD9emtlRR19ADSCUCDcjMTVDwcHTYJsVUszBpMJeu89YMYMa2ZFwj4ANwDshTlKZnw8EBb2J9DSjMDkEti82ZrZUMdqAC9DypqvL3D7dtWxQ0KAoqLKxj89XfyXxNmflRgVqlLvqKVmW3pnGv1Tsy4QaEE/v+0aZn/tZ0ubZtXV3g2WBiw13JCcTBM3xDAeNFfJkX5+PxBwZdu2s/T6305OttzmWrYGixCUlkYmJVmuBiUmai6FBMCmTZvWbsnYCAwRZFIbtHmzUK9LS80Tp+rthbrKe+zYMWzf7oYTJ06Y94C/IEwi6McfNcmRy01fBVcVKRIUCgXmzJmDsWPHwt3d3bRE/0YwydTTty+gXMAAV1dgqNY0s2itxQEREUL7SkwE/vgDyMurupP42Wef4cSJE5ikbV6uBwAjbpnr/WXXGqr0n13vN9vGUe9x0cZRT5CNo54gG0c9QTaOeoJsHPUE2TjqCbJx1BNk4/jbOjYHgIyMDABAVFQU5PKqv1VDzstNcWyek5ODK1euwM3NDffv34evry9KS0sRGhoKBweHesfmVQGQHJtv0hkALCsjb94UzsrVh0k6dya/+46cN093u+t27cglS8TQzH//Sx47Rl68SE6e/AFFUZPt28epnhseHi7l48/lVNbSUB93unOHPHGCTEsrUY5BndYoZDc30stLeLXy99cdx3J0JEePJt94Q31V3j0CoyiTtWVk5CTGxp5l9+5k27Zk48aks/McAp3o5KSgWM+0nmItLfjEE0/89QmKi7tIF5dpfPjhz7h3717u2rWbCxdu4eefZzI8vFRZiAUESmhnR0ZGkl26HKSPT2sdAry8yLw8sqKCXLBAOn+B0kIxqZaVlJBAFoGJBCIJPEsvr/2cNGkSAwIC+MUXX6jyl5OTw0ceaUdHxygKF2ekWL75F3VsLuH06dMMCUlSvuQAAiNoZ9eZcnl3Ojn1pptbGIWTc1/lrxsBT0ZHR9PPz499+87UO9JLkgcPipoRFvYy5XI3Sj6AevfuzeXLVxKYTsCPYWHv08MjnV27VhAgp04tZHp6Ops1a6ba07ukhExMVLBx4zUUPoV20NWVXLmy4K+xNYAkooYOJTMyyIULyRdfJP39xynJOUDJhYuHB3n7trgvO5sULstuKq/f5IABfzA9PZ07dpwnIESVNM/cx0fcl5srjp2cplOpKBG4y8ceO8K5c+cSiCMwkrt2HVflsbCwkO3bTyQAvvzyLNrZ2TEpKYl5ecXs109stlFURL755jKGhMRy48bK1d51RpCpfrOrQk4O2bWr5hfu6UmOGkV++SW5YcMFJUGf6NQAUrNWuLiI//fukeXlldckTJsmjouKKq/NmPEpfZWz7bX9eZ89q/mOLVt2UhNZq9m7dwrt7PwJgHFx+yltglxSUsKYmBjGx8er3MTUCUFpaZWF4upqmCSphsTHk2vWkK++SrZqJdoDbRGkvjohLY0EwNatV6lcN0tuCJ5+Whzn5FTGBzS3oVbzT0iFQvM527adZ7dug+jr24Rjx+rOQnJ0JDt0qHQZLZN1IQC6uXkQuKWMd5/AWrq4+PPLL79UzVYqKSnh1KlT6ePjU3dtkLbf7LFjdeM8eEA++qhmvKAgcs4c8sgRzS9du4ZIBWpnF8QVK1awtFQc9+1L/vab+P/555rP0yZc1x/3VQr3miMoHC29rRSPuvnQDDmUXGs6OrajTKbph65ly3OMiopi48aNOXToUO7atYukWJ6fnp5edzVIfXKinZ1wbH75Mvn116K2eHjoFpJ6DdmwQZyLixOiDSBnzRLXpClg48dP4tSpU0mKvol6Wikpuh+AdnB3J2NiSEfHS8qv+QkC/6aX12WS4oPQdiIr2if1WjWBABgcLDwFx8bO0IibnCw20UhPT+esWbMIgIGBgTx06BDJOhJxkudf7YmO9vbkM8+QK1cKdVZ7V6uBAyvTAEg/v8rj554T5xo0EL+enmR4eDjT0nbxyy8rSdQOQUH6Z8aGhIh0b9++TW/v9gQSKHZOEcQ8+eQDCvV8D4GTVToqDAkR/ukePHjAd999l88//3yVEy1zcnI4ZcoUhoaGMjg4mDt27Kh9ggytvnN0rIwnianQUKFBSXFIoakBQqNSh2Z6pyhU6Lt6n6VeG7/6SpNcJydxvry8nLGxsRw3bhwTEhSqDu0nn3zCykbfXfkc/e/k41NMADx9+jTHjh3L2bNnGy2joqIifvzxx2zSpEntE6Re2NohMLAy3vvvi3Pl5eI4J0ccSz14T0/Rh1iyRJCom15HAk144cJdkrofhtRe3bsnjuPjxfEbb4jjigry0qVLdHFxYUlJCUlyy5Yt3L9/P5s3f01Jjq+yFt2jv78Qd9o1adCg2wTA8ePH083NjRs2bDC5rHbu3Gk+QcuWLWOPHj2YmJjImJgYphnRm3v31k9O06bid9MmEQ8Qti11iH3e9IcBA9SP76m+8PPnz6vu196vjqw0yUhTvisqxPGbb4qGWqjCcbx7Vzg2b9bsJbVn3NMr2tT9Zt+5c4dhYWF0dnbmjBkzWFRUZDJBJ0+eNJ+ggQMHEgAjIyM5f/58nlV2Au7evavzwMREUSBubvp9fkZHV9YOqX04cICMjdVPiqdnZcFu315ZCwFSJkskAIMiZf16EXfvXs3zkhF0/PjxbNu2LdPT00mSmzZtoqvrySo/En3i0xQsWLCAzzzzDD/88EOuXLmS06dP56OPPkqrqdkFBQWqL1hyq2Vs5bf0UlVd79SJ7N5dv5iS7pPasC++EMcREW24Z88evXksKxNxWrXSvZadnUPgUwLgpUuXVPn39tbNl7096excdb6MIScnR1VWw4cPZ0JCAseOHUsAbNmypfXaoHPnznH58uW8fv06SeOerxo2rPqaOvSJKakm/PFH5TnhzgxMSEjg7NmzOXjwYDo4OHDjxo2sqKhgRIS4R+rFS8jPz+dDDz1EF5eRBKZQoVDofFz+/lJboyBwicDvDAzMUW5ueJdfffUVly5dypkzZzIzM5OHDh3ihAkTVESsXbuWCxYsYElJCRUKBVNSUhgXF8fg4GBOnDiRZ86cIQB+//331iOIJBctWsTt27erCtaYaIiJ0T0nqbuGoE+s7NnzgMBojhr1MV944QWuWLFCVUjffHOAwDW+/PLX7NatG8PDwxkREcHGjRsTAF977TVevSrS3blT+rguEviOwi/DBArfDZ4UTgObUSbzYufOneni4sKuXbsyMjKSQUFBDAoKYlRUFMeNG8ejR49K254RAA8dOsSJEyeyuLiYJJmRkcHJkyfTy8uLAGq3o/rkk4bJIStVaze3yvPKClglli0T8W7e1L0mbeIhQaFQcPDgwQSEfE9MTOS6det48uRJ/vbbbzx+/Dhzc3PVNvgrJ3CePj7xFBrbswTeYljYbHbs+AOF2UYYYgcNus0tW7bw3r17BvMrKR7SkEOuVl8hOzubKSkp3LdvnzIPNuDYXKoly5eLY6VWq1fEaQMQA1/6cPu2uP7BB5Xnunc/ofqCjx07pjpfUlLCvLw8koLIb7/9lg0bNqYYhphAoEhDTMvlwjxUk+WW/fr1Y8+ePfX6KJ03b55tODZX1+JOnaIyM5qEXL9Og42u1GGVhhD0YdiwStJPnSKB8ezb921+//339PHxYceOHRkQEKAibcCAZ9mwYTTl8kcIHDEqlmuyFla59ZmGA/OqUCsEBQdX/XLS/5s3xe+OHZr3SiOX+pSxqrQwdUhiU9jGrhLw5ZUrV1hWRi5fnssePX6lGBX9nWLj28UUm9dqb66hK371tX3GUFhYyFGjRrFFixacOXOm0fhWJ8iQePP1pcrSbOiLlKwPkugjhTUaIAsLjedBqMf3KEZWp+jkw8GBfOklcv/+yn6VZsdXU7SZU4Mee+wxJiQkcN++fSoLhSFYnSBD6nVCAjl3rnGxIfVZJBsZlbnr1En/MysqxNhRWBgptqKZRCCAwDMUo6hiRk1VMKXNrOlS/wkTJlAmk1Emk/HKlStG41udoJrsjerrK7SdzMxMlWlEzAMQ40HSCKfkuF2hILduJaOitNNaQEfHQIoNNLI0CtgQtDue+vJXUygUCq5atYoODg4EwH/84x/MzMzkrl27eOTIER1TkNUJqsr2BghDp0Khu6ddYiKZmZlJABw3bpwqLWlIQQr6+k2xseTPP4sdFJs0acKff/6dgK5FvCqYsh21JZxkVFRU8NixYxw3bhzDwsLYtWtXRkRE0M3NjS+++CJXrFhhdB9VswkyJirIykY8PFxXbFy5ckVl0ysqEgNo2ml06iQUC3X/FqWlpYyKiuKWLVs0npWfL/6//77+/Bqr7TXZiqa6yMnJ4YwZMxgVFcVp06ZZlyBD7Y+/v4izerU41tdeFhSQL79snGRtTJ06lf3792d6ejkB8vDhymvqarc2jO1OWZseGXfu3GmZXSANoao9U6Wh7Fu3dAv62jXdL7lZM/KXX3RHRRs21P/cmJgYpZlEl0Spxnbtqnl+1Sr9H5GlRZupUCgU/PTTT61LUPv2+gnSJsDTk+zVS1d0HT+umZ4kZozVoI4dO/K5574jIAjXhtTBlYyrYgaQCDJZ5eweudw8jc0SsBpBxvo/AQG65/v1q1r9/fprEUfq95w9K47nz9eN++67yQTe4COPGHrxSoKffVZ/DZdCXTppshpBhoa2q9OeVGaUbN688riwsJADB4rz6vbJAwcO0MUlgMBpjblt2jh8WNzbo4cgRDI9ubqSTZpo5qtNG8N5syasRpC+2f9SCA/XJTA4uPLemzdvMjMzU3W8YgVVbRZJSbvh5ctiNy0nJyGz09LS2LBhQwJpfOEF/flSH0/Stgo0by7EnfYsn79cDTJFvdY+5+1def+TTz7JsLAw7t+/n3FxcQTi6esr9OibN28SAEeOHEmFQsFffiGBO+zX7xWGh4czMnInAU21WyKlKpugPi3NVnzQWYUgQ+ItJKRyJs2aNSL+zz+L45UrxfGiRYv47bffStOOCAwiAB4/fpx+fn4EwKysLH7yySfs2bMn5XJPAiN49uwdAmLIOzFR5KM6otYWvTRahSBtjUw9HDhQqcWpo1s3ce7OnXJ6enpy27ZtSmvCUcpkBxkbG8vz58/z2rVrLC4u5tKlS9m4cWsOHPgDlyy5rfEMY6YadY2tZUvbdqFpcYKkmTvamlBiombBNWqkeZ/6BHUJu3eL4wsXtOMq2KRJWzo57TCJCGOhrsWYIViUIGOqtbYJX/urlTSrhQulzOnWtAMHDjIq6lF6eMRQDElXnxDtKV91qaUZg8UISk7WfXFjQZ/pRCJRmmLVt6+ww23fTr7ySint7ILo67uAcXFlOrszmhqkGT1/mxpkqltmqZ0x1Chrr8UBRM8+JoYcPnwF27ePYUWFUNHS0sjWrU17trZFwFa0NGOwCEHaS871BbJyuUdVjfK9e2KWjnZNlFTwMWPG8KuvvtK4x9Cz69pMYwnUSg2ytyfffVf8f/VVzQxUVAhlYNQoQURcHNmli25NS0tLo4+PD49rGejqHZubQFBamljyZ4qokcz8586JcZmmTUUtmD1b07Cp3uMfNWoUAfDgwYN6X+LPIq5qArMJUl9vakpYuFBMkm/QgBw/njx6VLPXrw1pot+ECROsVgi2DLMJ0l5vaiyEh5ObN4vZPKZg06ZNbNKkid5Jfn8HGCKo2n6zXVyAhx7SvK69WayDAzBwoPg1hlu3buGTTz7BO++8Y9Ch0d8VJpXIoEHAqlXA2LHA6tXAs89qXu/SRTe+qfjggw8QGRmJ0aNHm37T3wj1js1tA/WOzf+sqBf6No56gmwc9QTZOOoJsnHUE2TjqCfIxlFPkI2jniAbR833YjYRf0Xn6AXFBSgsKYSnkye8nb0NOyY3E1Yn6NatWzhy5K9hMbpTfAfT9kzDnP1zAAD3He7j84TPMThisL+1nml1gv7suH73OtLOpmHTmU3Yd2UffF18VdcelD3Aj7//aNXn17dBenAu7xxm/ToL3RZ1Q6sFrbDn0h6MjhqNqxOv4h/R/4AMMgCAq4Mr+jbva9W81NcgiEHLjOsZ2HRmEzae2YjbRbcxOHwwpvacil7NesHRzlEV92L+RSS0SkCgWyD6Nu+LQeHVGFupAf62BJUryrHn0h5sPL0Rm85ugou9C56MeBIpA1PQpVEXyGW6wqWkvATLji/D/tH70dy3ea3k829FkNRmbDyzEVuztqKZdzM8GfEkfhj5AyL9IyGTyQzen3Y2DW0D29YaOcDfgKDbRbexNWsrNp7ZiPSL6egY3BFDwodgeq/peMjrIeMJqCHlaArGRI2xUk714y9J0JU7V7DpzCZsOrsJh68exmOhj6nEl5+rX43SvJh/EcduHMPmyM0Wzq1h/CUIIonTt05j4+mN2HhmI7ILshEfFo83Or+Bvs37wtXB1exnfP/b93im7TNwtne2QI5Nx5+WoDd/eBObzmxCQ/eGyHuQh+LyYgyJGILZfWaje5PusJdb7tXKFeX497F/Y/sz2y2Wpqn4UxFEEoevHcYb29/AoauHAADZBdl4vt3zWDRokdFGvqb44fwPaOTZCG0D21olfUOweYJIIjMnE6syV2F15mo42Dkg70GeRpzD1w5bjRwASDmWgjHta1c5kGCzloSsvCxM+3ka2nzdBvGr4lGhqMD6YetxZuwZxDSO0YjbrVE3q+Xjxr0b+G/2f/FU66es9gxDsKkadPnOZazJXIPVJ1fj2t1rGNZqmN6OY7miHAPDBuJiwUU0826G3dm7kV+UDx8XH4vnaclvS5AQmQAPJw+Lp20K6pSg99Lfw/pT69HEuwnuld7D2VtnMTRyKGb3mY2eTXrCTm6nc0/egzz8cvkX/PHmH6pCm7BjAp7Z8Ay2jNii956agiRSjqVg+ZPLLZZmjTJhzVDVQuR3fnqHmApVGJY6jCXlxt1HLjyykElrkzTOlZaXMnZxLJN/suzalN0Xd7PNgjYq77xVAQYmv5sb6qwNWn9qvcbxqdxTGkbJqrAqcxVGtBmhcc7BzgFrE9di+YnlOumaA8lyYE0FxBjqjKBgj2CN4ydaPmH0nquFV/G/G/9DXMs4nWsBbgHYMGwDXtn2CjJzMs3OX35RPrZmbcXIh0eanZY5qDOC7pbexdNtn0argFYAgFm/zjJ6z5qTazAkYkiVvfkOwR3wed/P8eSaJ5FflG9W/pYfX464lnE1Ng1ZDNaSnVLQ1wZdzL9I/1n+LKsQu15k52cTU8G2C9oalPUdv+3IH8//aDAOSY7fPp79l/dneUW50bj6oFAo2HZBW+66sMuk+PirtUEbTm/A4PDBKnNME+8m+GHkDziRcwJv73xb7z3n8s7hyp0r6NWsl9H0Z/eZjeLyYry/+/0a5e/ItSO4X3YfsU1ja3S/JVEnBK0/vR4JkQka5/o274tpvaZh1r5ZSDuTpnPPqsxVGNZ6mEk2NklpWHliJdadWlft/KUcTcHoqNF6B+1qHdaqmlLQFnFXC6/S558+VarUjy56lJgKZt3K0hA54V+Ec9/lfSaJHAkZ1zLoP8ufJ24a3z9Bwt2Su/T+pzevFl41+R78lUTcxtMbMSBsQJUq9d4X9gIAwr4Mw/3S+wCA3278hpKKEnRt1LVaz2of1B5z+83FkNVDTFYa1p5ci+6Nu+tomXWFWido/en1GBox1GCc0imlAAD3T91BEqsyV2F46+E16o+MfHgkBoYNxNMbnkaFosJo/JSjKXix/YvVfo61UKsE5d7PRcb1DPRr0c9gPAc7B1ybeA0A4DbDDXP2zUHG9YwaP3d239koKS8xqjSczDmJ7IJsvf2sukKtEpR2Ns3kEc4gjyD0aNIDReVFUECBnRd2Iik1qUbPtZfbY03iGqNKw6Jji/B8u+ctOthnLmqVoA2nN+hob/qw/PhyyD6SYc+lPRrn0y+k1/jZAW4B2PDUBry67VW9lgZpStULUS/U+BnWQK0RVFBcgL2X92JAywF6r++7sg++M30h+0iGZzcKRwxdQjQdMPQO7W1WHgwpDZvObMLDgQ/X6pQqU1BrBG3N2orYprEa4yqXCi4helE0ZB/J8Oi/H0V+cT6m9ZqG8vfLwQ+J2X1mw8vJC77OvkhslYjUpFSz8yEpDSPWj9BQGlKO2ZZyoIK19HcpdOjQgcm7kukxw4PxK+NZWFzIUZtGaQw1PLvhWd4pvqPTv3h92+uc9vM0k/sjpqKsooy9FvfiOzvfIUleuH2B/rP8WVRm+hbP6oAV+0FWd2QRFB7EG0/f0DnfKbgTVieuRqhPqN77FFSg0eeNsPv/diPcP9zi+cq9n4tO33XCU62fwq6LuxDsEYzNI2o2500mk2WQ7GjhLAKohRHVguICjeNm3s1wYfwFo/f9evlX+Lv6W4UcQCgN47uMx8QfJwIATuaexOazm60+Gb66sHob5O3srXGc1No0VTn1VCqSWtVMrTYVv+f/rvpfXF5s9bU+NYHVCQrxCEFy92S0DmiNcL9w3Lh3A8bEqoIKrDu1zmQyawr1PlltrPWpEazVuElB3Vh6r+Qeo76J4sy9hvcW3ZO9x+jYkKWQdiaNY7eNZdqZtBqnASsqCbXaZXZzdMPmEZvRJaULIv0jMTB8oN54tSHeJAwKH2Rz7Y46at1Y2sizETYM24DRm0fjxM0TOtcVVGD96fVWF29/FtTJiFSXRl3wr/7/wqDVg5B7P1fjmrRQN8I/oi6yZnOosyHDp9s+jafbPI2ha4eitKJUdT71ZO2Jtz8D6nRMd1rvaQhwDcCrW18FSaG9nV5XT5Aa6pQguUyOpU8uRcb1DMw9MFcl3iIDIusyWzaFOh/4cHd0x+YRm9E1pStkMhkaeTSyyR59XcEGpq0Ajb0aY3yX8bh29xoOXTuEEetHYPPZ2l0LaquwCYIA4ErhFdX/2nCx8meBzRD0pzC71AGsPtzwN3GOXqVjcnNR79jcxmEzIq4e+lFPkI2jniAbRz1BNo56gmwc/w8z01Fe2H5dkAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 4 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, axs = plt.subplots(4, subplot_kw={'projection': proj, 'transform': transf})\n",
    "\n",
    "ax = axs[0]\n",
    "ax.coastlines()\n",
    "\n",
    "plot1(None, None, lons0, lats0, triang0, c0='r', c1='b', transf=transf, ax=ax)\n",
    "plot1(None, None, lons1, lats1, triang1, c0='c', c1='r', transf=transf, ax=ax)\n",
    "\n",
    "intersectedFalse = pystare.intersection(hull0, hull1, multi_resolution=False)\n",
    "intersectedTrue = pystare.intersection(hull0, hull1, multi_resolution=True)\n",
    "\n",
    "ax = axs[1]\n",
    "ax.coastlines()\n",
    "\n",
    "lonsF, latsF, intmatF = pystare.triangulate_indices(intersectedFalse)\n",
    "triangF = tri.Triangulation(lonsF, latsF, intmatF)\n",
    "plot1(None, None, lonsF, latsF, triangF, c0='r', c1='b', transf=transf, ax=ax)\n",
    "\n",
    "ax = axs[2]\n",
    "ax.coastlines()\n",
    "lonsT, latsT, intmatT = pystare.triangulate_indices(intersectedTrue)\n",
    "triangT = tri.Triangulation(lonsT, latsT, intmatT)\n",
    "plot1(None, None, lonsT, latsT, triangT, c0='r', c1='b', transf=transf, ax=ax)\n",
    "\n",
    "r0 = pystare.core.srange()\n",
    "r0.add_intervals(hull0)\n",
    "\n",
    "r1 = pystare.core.srange()\n",
    "r1.add_intervals(hull1)\n",
    "\n",
    "r01 = pystare.core.srange()\n",
    "r01.add_intersect(r0, r1, False)\n",
    "n01 = r01.get_size_as_values()\n",
    "\n",
    "intersected = numpy.zeros([n01], dtype=numpy.int64)\n",
    "r01.copy_values(intersected)\n",
    "# See examples/test_intersect_single_res.py\n",
    "\n",
    "r01.purge()\n",
    "n01 = r01.get_size_as_values()\n",
    "\n",
    "r01.reset()\n",
    "r01.add_intersect(r0, r1, True)\n",
    "n01 = r01.get_size_as_values()\n",
    "\n",
    "print('r01 n01: ', n01)\n",
    "\n",
    "intersected = numpy.zeros([n01], dtype=numpy.int64)\n",
    "r01.copy_values(intersected)\n",
    "\n",
    "###??? Would intervals be different?\n",
    "\n",
    "lonsRT, latsRT, intmatRT = pystare.triangulate_indices(intersected)\n",
    "triangRT = tri.Triangulation(lonsRT, latsRT, intmatRT)\n",
    "ax = axs[3]\n",
    "\n",
    "ax.coastlines()\n",
    "\n",
    "lonsRT_, latsRT_, intmatRT_ = pystare.triangulate_indices(intersected[51:55])\n",
    "triangRT_ = tri.Triangulation(lonsRT_, latsRT_, intmatRT_)\n",
    "\n",
    "plot1(None, None, lonsRT_, latsRT_, triangRT_, c0='g', c1='g', transf=transf, ax=ax)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}