2022-07-27_caustique.ipynb

{
 "cells": [
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "# Caustics\n",
    "\n",
    "[![DOI](https://zenodo.org/badge/273226625.svg)](https://zenodo.org/badge/latestdoi/273226625)\n",
    "\n",
    "Caustics ([wikipedia](https://en.wikipedia.org/wiki/Caustic_(optics))) are luminous patterns which are resulting from the superposition of smoothly deviated light rays. It is the heart-shaped pattern in your cup of coffee which is formed as the rays of from the sun are reflected on the cup's surface. It is also the wiggly patterns of light that you will see on the floor of a pool as the sun's light is *refracted* at the surface of the water. Here we will simulate this particular physical phenomenon. Simply because they are mesmerizingly beautiful, but also as it is of interest in visual neuroscience. Indeed, it speaks to how images are formed (more on this later), hence how the brain may understand images.\n",
    "\n",
    "In [this post](https://laurentperrinet.github.io/sciblog/posts/2020-06-19-caustic-optics.html), I have developed a simple formalism to generate such patterns, with the paradoxical result that it is *very* simple to code yet generates patterns with great complexity, such as:\n",
    "\n",
    "<BR>\n",
    "<center>\n",
    "<img src=\"caustique.gif\" width=\"100%\"/>\n",
    "</center>\n",
    "<BR>\n",
    "\n",
    "This is joint work with artist [Etienne Rey](https://laurentperrinet.github.io/authors/etienne-rey/), in which I especially follow the ideas put forward in the series [Turbulence](http://ondesparalleles.org/projets/turbulences/).\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "## installation\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:52.285561Z",
     "iopub.status.busy": "2022-07-28T16:28:52.284882Z",
     "iopub.status.idle": "2022-07-28T16:28:57.766264Z",
     "shell.execute_reply": "2022-07-28T16:28:57.766905Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Requirement already satisfied: pip in /usr/local/lib/python3.9/site-packages (from -r requirements.txt (line 4)) (22.2.1)\n",
      "Requirement already satisfied: matplotlib in /usr/local/lib/python3.9/site-packages (from -r requirements.txt (line 5)) (3.5.2)\n",
      "Requirement already satisfied: numpy in /usr/local/lib/python3.9/site-packages (from -r requirements.txt (line 6)) (1.23.1)\n",
      "Requirement already satisfied: MotionClouds in /usr/local/lib/python3.9/site-packages (from -r requirements.txt (line 7)) (20200212)\n",
      "Requirement already satisfied: imageio in /usr/local/lib/python3.9/site-packages (from -r requirements.txt (line 8)) (2.20.0)\n",
      "Requirement already satisfied: pygifsicle in /usr/local/lib/python3.9/site-packages (from -r requirements.txt (line 9)) (1.0.5)\n",
      "Requirement already satisfied: lambda2color in /usr/local/lib/python3.9/site-packages (from -r requirements.txt (line 10)) (1.0)\n",
      "Requirement already satisfied: pillow>=6.2.0 in /usr/local/lib/python3.9/site-packages (from matplotlib->-r requirements.txt (line 5)) (9.0.1)\n",
      "Requirement already satisfied: packaging>=20.0 in /usr/local/lib/python3.9/site-packages (from matplotlib->-r requirements.txt (line 5)) (20.9)\n",
      "Requirement already satisfied: python-dateutil>=2.7 in /usr/local/lib/python3.9/site-packages (from matplotlib->-r requirements.txt (line 5)) (2.8.2)\n",
      "Requirement already satisfied: pyparsing>=2.2.1 in /usr/local/lib/python3.9/site-packages (from matplotlib->-r requirements.txt (line 5)) (2.4.7)\n",
      "Requirement already satisfied: kiwisolver>=1.0.1 in /usr/local/lib/python3.9/site-packages (from matplotlib->-r requirements.txt (line 5)) (1.3.1)\n",
      "Requirement already satisfied: cycler>=0.10 in /usr/local/lib/python3.9/site-packages (from matplotlib->-r requirements.txt (line 5)) (0.10.0)\n",
      "Requirement already satisfied: fonttools>=4.22.0 in /usr/local/lib/python3.9/site-packages (from matplotlib->-r requirements.txt (line 5)) (4.30.0)\n",
      "Requirement already satisfied: six in /usr/local/lib/python3.9/site-packages (from cycler>=0.10->matplotlib->-r requirements.txt (line 5)) (1.16.0)\n",
      "Note: you may need to restart the kernel to use updated packages.\n"
     ]
    }
   ],
   "source": [
    "%pip install --upgrade -r requirements.txt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:57.771698Z",
     "iopub.status.busy": "2022-07-28T16:28:57.770961Z",
     "iopub.status.idle": "2022-07-28T16:28:57.772807Z",
     "shell.execute_reply": "2022-07-28T16:28:57.773187Z"
    }
   },
   "outputs": [],
   "source": [
    "figpath = '2022-07-30_caustique'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:57.777095Z",
     "iopub.status.busy": "2022-07-28T16:28:57.776595Z",
     "iopub.status.idle": "2022-07-28T16:28:58.434961Z",
     "shell.execute_reply": "2022-07-28T16:28:58.435353Z"
    }
   },
   "outputs": [],
   "source": [
    "from IPython.display import display\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# source code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:58.518062Z",
     "iopub.status.busy": "2022-07-28T16:28:58.517315Z",
     "iopub.status.idle": "2022-07-28T16:28:58.535632Z",
     "shell.execute_reply": "2022-07-28T16:28:58.535998Z"
    }
   },
   "outputs": [],
   "source": [
    "import os\n",
    "import datetime\n",
    "import numpy as np\n",
    "import MotionClouds as mc\n",
    "from lambda2color import Lambda2color, xyz_from_xy\n",
    "\n",
    "def init(args=[], figpath=figpath, PRECISION=7):\n",
    "\n",
    "    import argparse\n",
    "    if figpath is None:\n",
    "        import datetime\n",
    "        date = datetime.datetime.now().date().isoformat()\n",
    "        figpath = f'{date}_caustique'\n",
    "\n",
    "    parser = argparse.ArgumentParser()\n",
    "    parser.add_argument(\"--tag\", type=str, default='caustique', help=\"Tag\")\n",
    "    parser.add_argument(\"--figpath\", type=str, default=figpath, help=\"Folder to store images\")\n",
    "    parser.add_argument(\"--vext\", type=str, default='mp4', help=\"Folder to store images\")\n",
    "    parser.add_argument(\"--nx\", type=int, default=8*2**PRECISION, help=\"number of pixels (vertical)\")\n",
    "    parser.add_argument(\"--ny\", type=int, default=8*2**PRECISION, help=\"number of pixels (horizontal)\")\n",
    "    parser.add_argument(\"--nframe\", type=int, default=5*2**PRECISION, help=\"number of frames\")\n",
    "    parser.add_argument(\"--bin_dens\", type=int, default=2, help=\"relative bin density\")\n",
    "    parser.add_argument(\"--bin_spectrum\", type=int, default=6, help=\"bin spacing in spectrum\")\n",
    "    parser.add_argument(\"--seed\", type=int, default=42, help=\"seed for RNG\")\n",
    "    parser.add_argument(\"--H\", type=float, default=6., help=\"depth of the pool\")\n",
    "    parser.add_argument(\"--variation\", type=float, default=.4, help=\"variation of diffraction index: http://www.philiplaven.com/p20.html 1.40 at 400 nm and 1.37 at 700nm makes a 2% variation\")\n",
    "    parser.add_argument(\"--sf_0\", type=float, default=0.002, help=\"sf\")\n",
    "    parser.add_argument(\"--B_sf\", type=float, default=0.004, help=\"bandwidth in sf\")\n",
    "    parser.add_argument(\"--V_Y\", type=float, default=0.3, help=\"horizontal speed\")\n",
    "    parser.add_argument(\"--V_X\", type=float, default=0.3, help=\"vertical speed\")\n",
    "    parser.add_argument(\"--B_V\", type=float, default=4.0, help=\"bandwidth in speed\")\n",
    "    parser.add_argument(\"--theta\", type=float, default=2*np.pi*(2-1.61803), help=\"angle with the horizontal\")\n",
    "    parser.add_argument(\"--B_theta\", type=float, default=np.pi/3, help=\"bandwidth in theta\")\n",
    "    parser.add_argument(\"--min_lum\", type=float, default=.2, help=\"diffusion level for the rendering\")\n",
    "    parser.add_argument(\"--gamma\", type=float, default=2.8, help=\"Gamma exponant to convert luminosity to luminance\")\n",
    "    parser.add_argument(\"--fps\", type=float, default=18, help=\"frames per second\")\n",
    "    parser.add_argument(\"--multispectral\", type=bool, default=True, help=\"Compute caustics on the full spectrogram.\")\n",
    "    parser.add_argument(\"--cache\", type=bool, default=False, help=\"Cache intermediate output.\")\n",
    "    parser.add_argument(\"--verbose\", type=bool, default=False, help=\"Displays more verbose output.\")\n",
    "\n",
    "    opt = parser.parse_args(args=args)\n",
    "\n",
    "    if opt.verbose:\n",
    "        print(opt)\n",
    "    return opt\n",
    "\n",
    "def make_gif(gifname, fnames, fps, do_delete=True):\n",
    "    import imageio\n",
    "\n",
    "    with imageio.get_writer(gifname, mode='I', fps=fps) as writer:\n",
    "        for fname in fnames:\n",
    "            writer.append_data(imageio.imread(fname))\n",
    "\n",
    "    from pygifsicle import optimize\n",
    "    optimize(str(gifname))\n",
    "    if do_delete: \n",
    "        for fname in fnames: os.remove(fname)\n",
    "    return gifname\n",
    "\n",
    "# https://moviepy.readthedocs.io/en/latest/getting_started/videoclips.html#imagesequenceclip\n",
    "def make_mp4(mp4name, fnames, fps, do_delete=True):\n",
    "    import moviepy.editor as mpy\n",
    "    clip = mpy.ImageSequenceClip(fnames, fps=fps)\n",
    "    clip.write_videofile(mp4name, fps=fps, codec='libx264', verbose=False, logger=None)\n",
    "    if do_delete: \n",
    "        for fname in fnames: os.remove(fname)\n",
    "    return mp4name\n",
    "\n",
    "class Caustique:\n",
    "    def __init__(self, opt):\n",
    "        \"\"\"\n",
    "        Image coordinates follow 'ij' indexing, that is,\n",
    "        * their origin at the top left,\n",
    "        * the X axis is vertical and goes \"down\",\n",
    "        * the Y axis is horizontal and goes \"right\".\n",
    "\n",
    "        \"\"\"\n",
    "        self.mc = mc\n",
    "        self.ratio = opt.ny/opt.nx # ratio between height and width (>1 for portrait, <1 for landscape)\n",
    "        X = np.linspace(0, 1, opt.nx, endpoint=False) # vertical\n",
    "        Y = np.linspace(0, self.ratio, opt.ny, endpoint=False) # horizontal\n",
    "        self.xv, self.yv = np.meshgrid(X, Y, indexing='ij')\n",
    "        self.opt = opt\n",
    "        # https://stackoverflow.com/questions/16878315/what-is-the-right-way-to-treat-python-argparse-namespace-as-a-dictionary\n",
    "        self.d = vars(opt)\n",
    "        os.makedirs(self.opt.figpath, exist_ok=True)\n",
    "        self.cachepath = os.path.join('/tmp', self.opt.figpath)\n",
    "        if opt.verbose: print(f'{self.cachepath=}')\n",
    "        os.makedirs(self.cachepath, exist_ok=True)\n",
    "\n",
    "        # a standard white:\n",
    "        illuminant_D65 = xyz_from_xy(0.3127, 0.3291), \n",
    "        illuminant_sun = xyz_from_xy(0.325998, 0.335354)\n",
    "        # color conversion class\n",
    "        self.cs_srgb = Lambda2color(red=xyz_from_xy(0.64, 0.33),\n",
    "                               green=xyz_from_xy(0.30, 0.60),\n",
    "                               blue=xyz_from_xy(0.15, 0.06),\n",
    "                               white=illuminant_sun)\n",
    "    def wave(self):\n",
    "        filename = f'{self.cachepath}/{self.opt.tag}_wave.npy'\n",
    "        if os.path.isfile(filename):\n",
    "            z = np.load(filename)\n",
    "        else:\n",
    "            # A simplistic model of a wave using https://github.com/NeuralEnsemble/MotionClouds\n",
    "            import MotionClouds as mc\n",
    "            fx, fy, ft = mc.get_grids(self.opt.nx, self.opt.ny, self.opt.nframe)\n",
    "            env = mc.envelope_gabor(fx, fy, ft, V_X=self.opt.V_Y, V_Y=self.opt.V_X, B_V=self.opt.B_V,\n",
    "                                    sf_0=self.opt.sf_0, B_sf=self.opt.B_sf,\n",
    "                                    theta=self.opt.theta, B_theta=self.opt.B_theta)\n",
    "            z = mc.rectif(mc.random_cloud(env, seed=self.opt.seed))\n",
    "            if self.opt.cache: np.save(filename, z)\n",
    "        return z\n",
    "\n",
    "    def transform(self, z_, modulation=1.):\n",
    "        xv, yv = self.xv.copy(), self.yv.copy()\n",
    "\n",
    "        dzdx = z_ - np.roll(z_, 1, axis=0)\n",
    "        dzdy = z_ - np.roll(z_, 1, axis=1)\n",
    "        xv = xv + modulation * self.opt.H * dzdx\n",
    "        yv = yv + modulation * self.opt.H * dzdy\n",
    "\n",
    "        xv = np.mod(xv, 1)\n",
    "        yv = np.mod(yv, self.ratio)\n",
    "\n",
    "        return xv, yv\n",
    "\n",
    "    def do_raytracing(self, z):\n",
    "        filename = f'{self.cachepath}/{self.opt.tag}_hist.npy'\n",
    "        if os.path.isfile(filename):\n",
    "            hist = np.load(filename)\n",
    "        else:\n",
    "            binsx, binsy = self.opt.nx//self.opt.bin_dens, self.opt.ny//self.opt.bin_dens\n",
    "        \n",
    "            if self.opt.multispectral:\n",
    "                N_wavelengths = len(self.cs_srgb.cmf[:, 0])\n",
    "                \n",
    "                hist = np.zeros((binsx, binsy, self.opt.nframe, N_wavelengths//self.opt.bin_spectrum + 1))\n",
    "                for ii_wavelength, i_wavelength in enumerate(range(0, N_wavelengths, self.opt.bin_spectrum)):\n",
    "                    modulation = 1. + self.opt.variation/2 - self.opt.variation*i_wavelength/N_wavelengths\n",
    "                    # print(i_wavelength, N_wavelengths, modulation)\n",
    "                    for i_frame in range(self.opt.nframe):\n",
    "                        xv, yv = self.transform(z[:, :, i_frame], modulation=modulation)\n",
    "                        hist_, edge_x, edge_y = np.histogram2d(xv.ravel(), yv.ravel(),\n",
    "                                                               bins=[binsx, binsy],\n",
    "                                                               range=[[0, 1], [0, self.ratio]],\n",
    "                                                               density=True)\n",
    "                        hist[:, :, i_frame, ii_wavelength] = hist_\n",
    "                hist /= hist.max()\n",
    "            else:\n",
    "                hist = np.zeros((binsx, binsy, self.opt.nframe))\n",
    "                for i_frame in range(self.opt.nframe):\n",
    "                    xv, yv = self.transform(z[:, :, i_frame])\n",
    "                    hist_, edge_x, edge_y = np.histogram2d(xv.ravel(), yv.ravel(),\n",
    "                                                           bins=[binsx, binsy],\n",
    "                                                           range=[[0, 1], [0, self.ratio]],\n",
    "                                                           density=True)\n",
    "                #hist /= hist.max()\n",
    "            if self.opt.cache: np.save(filename, hist)\n",
    "        return hist\n",
    "    \n",
    "\n",
    "    def plot(self, z, do_color=True, videoname=None, dpi=50):\n",
    "        hist = self.do_raytracing(z)\n",
    "\n",
    "        if videoname is None:\n",
    "            videoname=f'{self.opt.figpath}/{self.opt.tag}.{self.opt.vext}'\n",
    "\n",
    "        subplotpars = matplotlib.figure.SubplotParams(left=0., right=1., bottom=0., top=1., wspace=0., hspace=0.,)\n",
    "\n",
    "        if self.opt.multispectral:\n",
    "            N_wavelengths = len(self.cs_srgb.cmf[:, 0])\n",
    "            # multiply by the spectrum of the sky\n",
    "            if True:\n",
    "                wavelengths = self.cs_srgb.cmf[::self.opt.bin_spectrum, 0]*1e-9\n",
    "                intensity5800 = planck(wavelengths, 5800.)\n",
    "                scatter = scattering(wavelengths)\n",
    "                spectrum_sky = intensity5800 * scatter\n",
    "                hist = hist * spectrum_sky[None, None, None, :]\n",
    "                hist /= hist.max()\n",
    "\n",
    "            if False:\n",
    "                # some magic to only get the hue\n",
    "                hist_max = hist.max(axis=-1)[:, :, :, None]*np.ones((1, 1, 1, N_wavelengths//self.opt.bin_spectrum + 1))\n",
    "                #print(hist.shape, hist.max(axis=-1).shape, hist_max.shape)\n",
    "                hist[hist_max==0] = 1\n",
    "                hist_max[hist_max==0] = 1\n",
    "                hist /= hist_max\n",
    "                #hist -= hist.min()\n",
    "            else:\n",
    "                hist /= hist.max()\n",
    "            \n",
    "            #image_rgb = self.cs_srgb.spec_to_rgb(hist)\n",
    "            image_rgb = np.zeros((self.opt.nx//self.opt.bin_dens,  self.opt.ny//self.opt.bin_dens, 3, self.opt.nframe))\n",
    "            for ii_wavelength, i_wavelength in enumerate(range(0, N_wavelengths, self.opt.bin_spectrum)):\n",
    "                spec = np.zeros((N_wavelengths))\n",
    "                spec[i_wavelength] = 1\n",
    "                rgb = self.cs_srgb.spec_to_rgb(spec)\n",
    "                # print(hist.shape, image_rgb.shape, rgb.shape)\n",
    "                image_rgb += hist[:, :, None, :, ii_wavelength] * rgb[None, None, :, None]\n",
    "            \n",
    "            # image_rgb -= image_rgb.min()\n",
    "            image_rgb /= image_rgb.max()\n",
    "\n",
    "\n",
    "        fnames = []\n",
    "        for i_frame in range(self.opt.nframe):\n",
    "            fig, ax = plt.subplots(figsize=(self.opt.ny/self.opt.bin_dens/dpi, self.opt.nx/self.opt.bin_dens/dpi), subplotpars=subplotpars)\n",
    "            if self.opt.multispectral:\n",
    "                ax.imshow(image_rgb[:, :, :, i_frame] ** (1/self.opt.gamma), vmin=0, vmax=1)\n",
    "            else:\n",
    "                if do_color:\n",
    "                    bluesky = np.array([0.268375, 0.283377]) # xyz\n",
    "                    sun = np.array([0.325998, 0.335354]) # xyz\n",
    "                    # ax.pcolormesh(edge_y, edge_x, hist[:, :, i_frame], vmin=0, vmax=1, cmap=plt.cm.Blues_r)\n",
    "                    # https://en.wikipedia.org/wiki/CIE_1931_color_space#Mixing_colors_specified_with_the_CIE_xy_chromaticity_diagram\n",
    "                    L1 = 1 - hist[:, :, i_frame]\n",
    "                    L2 = hist[:, :, i_frame]\n",
    "                    image_denom = L1 / bluesky[1] + L2 / sun[1]\n",
    "                    image_x = (L1 * bluesky[0] / bluesky[1] + L2 * sun[0] / sun[1]) / image_denom\n",
    "                    image_y = (L1 + L2) / image_denom \n",
    "                    image_xyz = np.dstack((image_x, image_y, 1 - image_x - image_y))\n",
    "                    image_rgb = self.cs_srgb.xyz_to_rgb(image_xyz)\n",
    "                    image_L = self.opt.min_lum + (1-self.opt.min_lum)* L2 ** .61803\n",
    "                    ax.imshow(image_L[:, :, None]*image_rgb, vmin=0, vmax=1)\n",
    "\n",
    "                else:\n",
    "                    ax.imshow(1-image_L, vmin=0, vmax=1)\n",
    "\n",
    "            fname = f'{self.cachepath}/{self.opt.tag}_frame_{i_frame:04d}.png'\n",
    "            fig.savefig(fname, dpi=dpi)\n",
    "            fnames.append(fname)\n",
    "            plt.close('all')\n",
    "        if self.opt.vext == 'gif':\n",
    "            return make_gif(videoname, fnames, fps=self.opt.fps)\n",
    "        else:\n",
    "            return make_mp4(videoname, fnames, fps=self.opt.fps)\n",
    "\n",
    "    def show(self, videoname, width=1024):\n",
    "        if self.opt.vext == 'gif':\n",
    "            from IPython.display import Image, display\n",
    "            return display(Image(url=videoname, width=width))\n",
    "        else:\n",
    "            #import moviepy.editor as mpy\n",
    "            #return mpy.ipython_display(videoname, width=width)\n",
    "            # https://github.com/NeuralEnsemble/MotionClouds/blob/master/MotionClouds/MotionClouds.py#L858\n",
    "            from IPython.display import HTML, display\n",
    "            opts = ' loop=\"1\" autoplay=\"1\" controls '\n",
    "            html = HTML(f'<video {opts} src=\"{videoname}\" type=\"video/{self.opt.vext}\" width={width}\\>')\n",
    "            html.reload()\n",
    "            return display(html)\n",
    "    \n",
    "# borrowed from https://github.com/gummiks/gummiks.github.io/blob/master/scripts/astro/planck.py\n",
    "def planck(wav, T):\n",
    "    import scipy.constants as const\n",
    "    c = const.c # c = 3.0e+8\n",
    "    h = const.h # h = 6.626e-34\n",
    "    k = const.k # k = 1.38e-23\n",
    "    a = 2.0*h*c**2\n",
    "    b = h*c/(wav*k*T)\n",
    "    intensity = a / ( (wav**5) * (np.exp(b) - 1.0) )\n",
    "    return intensity\n",
    "\n",
    "def scattering(wav, a=0.005, p=1.3, b=0.45):\n",
    "    \"\"\"\n",
    "    b is  proportionate  to  the  column  density  of  aerosols  \n",
    "    along  the  path  of  sunlight,  from  outside  the  atmosphere \n",
    "    to  the  point  of  observation\n",
    "    \n",
    "    N_O3  is  the  ozone  column  density  along  the  path  of  sunlight,  \n",
    "    sigma_O3 is  the  wavelength dependent ozone absorption cross-section.\n",
    "    \n",
    "    \"\"\"\n",
    "    # converting wav in µm:\n",
    "    intensity = np.exp(-a/((wav/1e-6)**4)) # Rayleigh extinction by nitrogen\n",
    "    intensity *= (wav/1e-6)**-4\n",
    "    intensity *= np.exp(-b/((wav/1e-6)**p)) # Aerosols\n",
    "    return intensity"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "## parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:58.539730Z",
     "iopub.status.busy": "2022-07-28T16:28:58.539143Z",
     "iopub.status.idle": "2022-07-28T16:28:58.541363Z",
     "shell.execute_reply": "2022-07-28T16:28:58.542023Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Saving our simulations in=None\n"
     ]
    }
   ],
   "source": [
    "print(f'Saving our simulations in={figpath}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:58.545418Z",
     "iopub.status.busy": "2022-07-28T16:28:58.544922Z",
     "iopub.status.idle": "2022-07-28T16:28:58.547859Z",
     "shell.execute_reply": "2022-07-28T16:28:58.548435Z"
    }
   },
   "outputs": [],
   "source": [
    "opt = init()\n",
    "opt.figpath = figpath"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:58.555813Z",
     "iopub.status.busy": "2022-07-28T16:28:58.555174Z",
     "iopub.status.idle": "2022-07-28T16:28:58.557606Z",
     "shell.execute_reply": "2022-07-28T16:28:58.558030Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Namespace(tag='caustique', figpath=None, vext='mp4', nx=1024, ny=1024, nframe=640, bin_dens=2, bin_spectrum=6, seed=42, H=10.0, variation=0.4, sf_0=0.004, B_sf=0.004, V_Y=0.3, V_X=0.3, B_V=4.0, theta=2.399988291783386, B_theta=1.0471975511965976, min_lum=0.2, gamma=2.4, fps=18, multispectral=True, cache=False, verbose=False)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "opt"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "%rm -fr {figpath}"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "%rm -fr {figpath}/*.{opt.vext}"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "%rm -fr /tmp/{figpath} "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# une simple caustique"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "## a wave pattern"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:58.562728Z",
     "iopub.status.busy": "2022-07-28T16:28:58.562155Z",
     "iopub.status.idle": "2022-07-28T16:28:59.270442Z",
     "shell.execute_reply": "2022-07-28T16:28:59.270945Z"
    }
   },
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "expected str, bytes or os.PathLike object, not NoneType",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Input \u001b[0;32mIn [8]\u001b[0m, in \u001b[0;36m<cell line: 1>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0m c \u001b[38;5;241m=\u001b[39m \u001b[43mCaustique\u001b[49m\u001b[43m(\u001b[49m\u001b[43mopt\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m      2\u001b[0m wave_videoname \u001b[38;5;241m=\u001b[39m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mopt\u001b[38;5;241m.\u001b[39mfigpath\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m/\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mopt\u001b[38;5;241m.\u001b[39mtag\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m_wave\u001b[39m\u001b[38;5;124m'\u001b[39m\n\u001b[1;32m      3\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m os\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39misfile(\u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mwave_videoname\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;132;01m{\u001b[39;00mc\u001b[38;5;241m.\u001b[39mmc\u001b[38;5;241m.\u001b[39mvext\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m'\u001b[39m):\n",
      "Input \u001b[0;32mIn [4]\u001b[0m, in \u001b[0;36mCaustique.__init__\u001b[0;34m(self, opt)\u001b[0m\n\u001b[1;32m     86\u001b[0m \u001b[38;5;66;03m# https://stackoverflow.com/questions/16878315/what-is-the-right-way-to-treat-python-argparse-namespace-as-a-dictionary\u001b[39;00m\n\u001b[1;32m     87\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39md \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mvars\u001b[39m(opt)\n\u001b[0;32m---> 88\u001b[0m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mmakedirs\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mopt\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfigpath\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mexist_ok\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;28;43;01mTrue\u001b[39;49;00m\u001b[43m)\u001b[49m\n\u001b[1;32m     89\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath \u001b[38;5;241m=\u001b[39m os\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39mjoin(\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m/tmp\u001b[39m\u001b[38;5;124m'\u001b[39m, \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mopt\u001b[38;5;241m.\u001b[39mfigpath)\n\u001b[1;32m     90\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m opt\u001b[38;5;241m.\u001b[39mverbose: \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00m\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath\u001b[38;5;132;01m=}\u001b[39;00m\u001b[38;5;124m'\u001b[39m)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/os.py:210\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n\u001b[1;32m    200\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mmakedirs\u001b[39m(name, mode\u001b[38;5;241m=\u001b[39m\u001b[38;5;241m0o777\u001b[39m, exist_ok\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m):\n\u001b[1;32m    201\u001b[0m     \u001b[38;5;124;03m\"\"\"makedirs(name [, mode=0o777][, exist_ok=False])\u001b[39;00m\n\u001b[1;32m    202\u001b[0m \n\u001b[1;32m    203\u001b[0m \u001b[38;5;124;03m    Super-mkdir; create a leaf directory and all intermediate ones.  Works like\u001b[39;00m\n\u001b[0;32m   (...)\u001b[0m\n\u001b[1;32m    208\u001b[0m \n\u001b[1;32m    209\u001b[0m \u001b[38;5;124;03m    \"\"\"\u001b[39;00m\n\u001b[0;32m--> 210\u001b[0m     head, tail \u001b[38;5;241m=\u001b[39m \u001b[43mpath\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msplit\u001b[49m\u001b[43m(\u001b[49m\u001b[43mname\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    211\u001b[0m     \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m tail:\n\u001b[1;32m    212\u001b[0m         head, tail \u001b[38;5;241m=\u001b[39m path\u001b[38;5;241m.\u001b[39msplit(head)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/posixpath.py:103\u001b[0m, in \u001b[0;36msplit\u001b[0;34m(p)\u001b[0m\n\u001b[1;32m    100\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21msplit\u001b[39m(p):\n\u001b[1;32m    101\u001b[0m     \u001b[38;5;124;03m\"\"\"Split a pathname.  Returns tuple \"(head, tail)\" where \"tail\" is\u001b[39;00m\n\u001b[1;32m    102\u001b[0m \u001b[38;5;124;03m    everything after the final slash.  Either part may be empty.\"\"\"\u001b[39;00m\n\u001b[0;32m--> 103\u001b[0m     p \u001b[38;5;241m=\u001b[39m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfspath\u001b[49m\u001b[43m(\u001b[49m\u001b[43mp\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    104\u001b[0m     sep \u001b[38;5;241m=\u001b[39m _get_sep(p)\n\u001b[1;32m    105\u001b[0m     i \u001b[38;5;241m=\u001b[39m p\u001b[38;5;241m.\u001b[39mrfind(sep) \u001b[38;5;241m+\u001b[39m \u001b[38;5;241m1\u001b[39m\n",
      "\u001b[0;31mTypeError\u001b[0m: expected str, bytes or os.PathLike object, not NoneType"
     ]
    }
   ],
   "source": [
    "c = Caustique(opt)\n",
    "wave_videoname = f'{opt.figpath}/{opt.tag}_wave'\n",
    "if not os.path.isfile(f'{wave_videoname}{c.mc.vext}'):\n",
    "    print(f'Doing {wave_videoname}{c.mc.vext}')\n",
    "    z = c.wave()\n",
    "    c.mc.anim_save(z.swapaxes(0, 1), wave_videoname)\n",
    "c.show(f'{wave_videoname}{c.mc.vext}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "## generating caustics"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:59.274870Z",
     "iopub.status.busy": "2022-07-28T16:28:59.274323Z",
     "iopub.status.idle": "2022-07-28T16:28:59.311684Z",
     "shell.execute_reply": "2022-07-28T16:28:59.312262Z"
    }
   },
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "expected str, bytes or os.PathLike object, not NoneType",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Input \u001b[0;32mIn [9]\u001b[0m, in \u001b[0;36m<cell line: 1>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0m c \u001b[38;5;241m=\u001b[39m \u001b[43mCaustique\u001b[49m\u001b[43m(\u001b[49m\u001b[43mopt\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m      2\u001b[0m z \u001b[38;5;241m=\u001b[39m c\u001b[38;5;241m.\u001b[39mwave()\n\u001b[1;32m      3\u001b[0m z\u001b[38;5;241m.\u001b[39mshape\n",
      "Input \u001b[0;32mIn [4]\u001b[0m, in \u001b[0;36mCaustique.__init__\u001b[0;34m(self, opt)\u001b[0m\n\u001b[1;32m     86\u001b[0m \u001b[38;5;66;03m# https://stackoverflow.com/questions/16878315/what-is-the-right-way-to-treat-python-argparse-namespace-as-a-dictionary\u001b[39;00m\n\u001b[1;32m     87\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39md \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mvars\u001b[39m(opt)\n\u001b[0;32m---> 88\u001b[0m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mmakedirs\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mopt\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfigpath\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mexist_ok\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;28;43;01mTrue\u001b[39;49;00m\u001b[43m)\u001b[49m\n\u001b[1;32m     89\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath \u001b[38;5;241m=\u001b[39m os\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39mjoin(\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m/tmp\u001b[39m\u001b[38;5;124m'\u001b[39m, \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mopt\u001b[38;5;241m.\u001b[39mfigpath)\n\u001b[1;32m     90\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m opt\u001b[38;5;241m.\u001b[39mverbose: \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00m\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath\u001b[38;5;132;01m=}\u001b[39;00m\u001b[38;5;124m'\u001b[39m)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/os.py:210\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n\u001b[1;32m    200\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mmakedirs\u001b[39m(name, mode\u001b[38;5;241m=\u001b[39m\u001b[38;5;241m0o777\u001b[39m, exist_ok\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m):\n\u001b[1;32m    201\u001b[0m     \u001b[38;5;124;03m\"\"\"makedirs(name [, mode=0o777][, exist_ok=False])\u001b[39;00m\n\u001b[1;32m    202\u001b[0m \n\u001b[1;32m    203\u001b[0m \u001b[38;5;124;03m    Super-mkdir; create a leaf directory and all intermediate ones.  Works like\u001b[39;00m\n\u001b[0;32m   (...)\u001b[0m\n\u001b[1;32m    208\u001b[0m \n\u001b[1;32m    209\u001b[0m \u001b[38;5;124;03m    \"\"\"\u001b[39;00m\n\u001b[0;32m--> 210\u001b[0m     head, tail \u001b[38;5;241m=\u001b[39m \u001b[43mpath\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msplit\u001b[49m\u001b[43m(\u001b[49m\u001b[43mname\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    211\u001b[0m     \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m tail:\n\u001b[1;32m    212\u001b[0m         head, tail \u001b[38;5;241m=\u001b[39m path\u001b[38;5;241m.\u001b[39msplit(head)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/posixpath.py:103\u001b[0m, in \u001b[0;36msplit\u001b[0;34m(p)\u001b[0m\n\u001b[1;32m    100\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21msplit\u001b[39m(p):\n\u001b[1;32m    101\u001b[0m     \u001b[38;5;124;03m\"\"\"Split a pathname.  Returns tuple \"(head, tail)\" where \"tail\" is\u001b[39;00m\n\u001b[1;32m    102\u001b[0m \u001b[38;5;124;03m    everything after the final slash.  Either part may be empty.\"\"\"\u001b[39;00m\n\u001b[0;32m--> 103\u001b[0m     p \u001b[38;5;241m=\u001b[39m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfspath\u001b[49m\u001b[43m(\u001b[49m\u001b[43mp\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    104\u001b[0m     sep \u001b[38;5;241m=\u001b[39m _get_sep(p)\n\u001b[1;32m    105\u001b[0m     i \u001b[38;5;241m=\u001b[39m p\u001b[38;5;241m.\u001b[39mrfind(sep) \u001b[38;5;241m+\u001b[39m \u001b[38;5;241m1\u001b[39m\n",
      "\u001b[0;31mTypeError\u001b[0m: expected str, bytes or os.PathLike object, not NoneType"
     ]
    }
   ],
   "source": [
    "c = Caustique(opt)\n",
    "z = c.wave()\n",
    "z.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:59.316187Z",
     "iopub.status.busy": "2022-07-28T16:28:59.315685Z",
     "iopub.status.idle": "2022-07-28T16:28:59.354474Z",
     "shell.execute_reply": "2022-07-28T16:28:59.354880Z"
    }
   },
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "expected str, bytes or os.PathLike object, not NoneType",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Input \u001b[0;32mIn [10]\u001b[0m, in \u001b[0;36m<cell line: 2>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m videoname \u001b[38;5;241m=\u001b[39m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mopt\u001b[38;5;241m.\u001b[39mfigpath\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m/\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mopt\u001b[38;5;241m.\u001b[39mtag\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m.\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mopt\u001b[38;5;241m.\u001b[39mvext\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m'\u001b[39m\n\u001b[1;32m      2\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m os\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39misfile(videoname):\n\u001b[0;32m----> 3\u001b[0m     c \u001b[38;5;241m=\u001b[39m \u001b[43mCaustique\u001b[49m\u001b[43m(\u001b[49m\u001b[43mopt\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m      4\u001b[0m     z \u001b[38;5;241m=\u001b[39m c\u001b[38;5;241m.\u001b[39mwave()\n\u001b[1;32m      5\u001b[0m     videoname \u001b[38;5;241m=\u001b[39m c\u001b[38;5;241m.\u001b[39mplot(z)\n",
      "Input \u001b[0;32mIn [4]\u001b[0m, in \u001b[0;36mCaustique.__init__\u001b[0;34m(self, opt)\u001b[0m\n\u001b[1;32m     86\u001b[0m \u001b[38;5;66;03m# https://stackoverflow.com/questions/16878315/what-is-the-right-way-to-treat-python-argparse-namespace-as-a-dictionary\u001b[39;00m\n\u001b[1;32m     87\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39md \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mvars\u001b[39m(opt)\n\u001b[0;32m---> 88\u001b[0m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mmakedirs\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mopt\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfigpath\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mexist_ok\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;28;43;01mTrue\u001b[39;49;00m\u001b[43m)\u001b[49m\n\u001b[1;32m     89\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath \u001b[38;5;241m=\u001b[39m os\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39mjoin(\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m/tmp\u001b[39m\u001b[38;5;124m'\u001b[39m, \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mopt\u001b[38;5;241m.\u001b[39mfigpath)\n\u001b[1;32m     90\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m opt\u001b[38;5;241m.\u001b[39mverbose: \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00m\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath\u001b[38;5;132;01m=}\u001b[39;00m\u001b[38;5;124m'\u001b[39m)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/os.py:210\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n\u001b[1;32m    200\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mmakedirs\u001b[39m(name, mode\u001b[38;5;241m=\u001b[39m\u001b[38;5;241m0o777\u001b[39m, exist_ok\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m):\n\u001b[1;32m    201\u001b[0m     \u001b[38;5;124;03m\"\"\"makedirs(name [, mode=0o777][, exist_ok=False])\u001b[39;00m\n\u001b[1;32m    202\u001b[0m \n\u001b[1;32m    203\u001b[0m \u001b[38;5;124;03m    Super-mkdir; create a leaf directory and all intermediate ones.  Works like\u001b[39;00m\n\u001b[0;32m   (...)\u001b[0m\n\u001b[1;32m    208\u001b[0m \n\u001b[1;32m    209\u001b[0m \u001b[38;5;124;03m    \"\"\"\u001b[39;00m\n\u001b[0;32m--> 210\u001b[0m     head, tail \u001b[38;5;241m=\u001b[39m \u001b[43mpath\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msplit\u001b[49m\u001b[43m(\u001b[49m\u001b[43mname\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    211\u001b[0m     \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m tail:\n\u001b[1;32m    212\u001b[0m         head, tail \u001b[38;5;241m=\u001b[39m path\u001b[38;5;241m.\u001b[39msplit(head)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/posixpath.py:103\u001b[0m, in \u001b[0;36msplit\u001b[0;34m(p)\u001b[0m\n\u001b[1;32m    100\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21msplit\u001b[39m(p):\n\u001b[1;32m    101\u001b[0m     \u001b[38;5;124;03m\"\"\"Split a pathname.  Returns tuple \"(head, tail)\" where \"tail\" is\u001b[39;00m\n\u001b[1;32m    102\u001b[0m \u001b[38;5;124;03m    everything after the final slash.  Either part may be empty.\"\"\"\u001b[39;00m\n\u001b[0;32m--> 103\u001b[0m     p \u001b[38;5;241m=\u001b[39m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfspath\u001b[49m\u001b[43m(\u001b[49m\u001b[43mp\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    104\u001b[0m     sep \u001b[38;5;241m=\u001b[39m _get_sep(p)\n\u001b[1;32m    105\u001b[0m     i \u001b[38;5;241m=\u001b[39m p\u001b[38;5;241m.\u001b[39mrfind(sep) \u001b[38;5;241m+\u001b[39m \u001b[38;5;241m1\u001b[39m\n",
      "\u001b[0;31mTypeError\u001b[0m: expected str, bytes or os.PathLike object, not NoneType"
     ]
    }
   ],
   "source": [
    "videoname = f'{opt.figpath}/{opt.tag}.{opt.vext}'\n",
    "if not os.path.isfile(videoname):\n",
    "    c = Caustique(opt)\n",
    "    z = c.wave()\n",
    "    videoname = c.plot(z)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:59.368027Z",
     "iopub.status.busy": "2022-07-28T16:28:59.367438Z",
     "iopub.status.idle": "2022-07-28T16:28:59.369837Z",
     "shell.execute_reply": "2022-07-28T16:28:59.370268Z"
    }
   },
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'c' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "Input \u001b[0;32mIn [11]\u001b[0m, in \u001b[0;36m<cell line: 1>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0m \u001b[43mc\u001b[49m\u001b[38;5;241m.\u001b[39mshow(videoname)\n",
      "\u001b[0;31mNameError\u001b[0m: name 'c' is not defined"
     ]
    }
   ],
   "source": [
    "c.show(videoname)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "## exploring parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:59.373680Z",
     "iopub.status.busy": "2022-07-28T16:28:59.373010Z",
     "iopub.status.idle": "2022-07-28T16:28:59.374895Z",
     "shell.execute_reply": "2022-07-28T16:28:59.375315Z"
    }
   },
   "outputs": [],
   "source": [
    "N_scan = 9\n",
    "base = 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:59.380520Z",
     "iopub.status.busy": "2022-07-28T16:28:59.379958Z",
     "iopub.status.idle": "2022-07-28T16:28:59.425627Z",
     "shell.execute_reply": "2022-07-28T16:28:59.426066Z"
    }
   },
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "expected str, bytes or os.PathLike object, not NoneType",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Input \u001b[0;32mIn [13]\u001b[0m, in \u001b[0;36m<cell line: 4>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m opt \u001b[38;5;241m=\u001b[39m init()\n\u001b[1;32m      2\u001b[0m opt\u001b[38;5;241m.\u001b[39mfigpath \u001b[38;5;241m=\u001b[39m figpath\n\u001b[0;32m----> 4\u001b[0m c \u001b[38;5;241m=\u001b[39m \u001b[43mCaustique\u001b[49m\u001b[43m(\u001b[49m\u001b[43mopt\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m      5\u001b[0m \u001b[38;5;66;03m# compute just once\u001b[39;00m\n\u001b[1;32m      6\u001b[0m z \u001b[38;5;241m=\u001b[39m c\u001b[38;5;241m.\u001b[39mwave()\n",
      "Input \u001b[0;32mIn [4]\u001b[0m, in \u001b[0;36mCaustique.__init__\u001b[0;34m(self, opt)\u001b[0m\n\u001b[1;32m     86\u001b[0m \u001b[38;5;66;03m# https://stackoverflow.com/questions/16878315/what-is-the-right-way-to-treat-python-argparse-namespace-as-a-dictionary\u001b[39;00m\n\u001b[1;32m     87\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39md \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mvars\u001b[39m(opt)\n\u001b[0;32m---> 88\u001b[0m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mmakedirs\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mopt\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfigpath\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mexist_ok\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;28;43;01mTrue\u001b[39;49;00m\u001b[43m)\u001b[49m\n\u001b[1;32m     89\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath \u001b[38;5;241m=\u001b[39m os\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39mjoin(\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m/tmp\u001b[39m\u001b[38;5;124m'\u001b[39m, \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mopt\u001b[38;5;241m.\u001b[39mfigpath)\n\u001b[1;32m     90\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m opt\u001b[38;5;241m.\u001b[39mverbose: \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00m\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath\u001b[38;5;132;01m=}\u001b[39;00m\u001b[38;5;124m'\u001b[39m)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/os.py:210\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n\u001b[1;32m    200\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mmakedirs\u001b[39m(name, mode\u001b[38;5;241m=\u001b[39m\u001b[38;5;241m0o777\u001b[39m, exist_ok\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m):\n\u001b[1;32m    201\u001b[0m     \u001b[38;5;124;03m\"\"\"makedirs(name [, mode=0o777][, exist_ok=False])\u001b[39;00m\n\u001b[1;32m    202\u001b[0m \n\u001b[1;32m    203\u001b[0m \u001b[38;5;124;03m    Super-mkdir; create a leaf directory and all intermediate ones.  Works like\u001b[39;00m\n\u001b[0;32m   (...)\u001b[0m\n\u001b[1;32m    208\u001b[0m \n\u001b[1;32m    209\u001b[0m \u001b[38;5;124;03m    \"\"\"\u001b[39;00m\n\u001b[0;32m--> 210\u001b[0m     head, tail \u001b[38;5;241m=\u001b[39m \u001b[43mpath\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msplit\u001b[49m\u001b[43m(\u001b[49m\u001b[43mname\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    211\u001b[0m     \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m tail:\n\u001b[1;32m    212\u001b[0m         head, tail \u001b[38;5;241m=\u001b[39m path\u001b[38;5;241m.\u001b[39msplit(head)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/posixpath.py:103\u001b[0m, in \u001b[0;36msplit\u001b[0;34m(p)\u001b[0m\n\u001b[1;32m    100\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21msplit\u001b[39m(p):\n\u001b[1;32m    101\u001b[0m     \u001b[38;5;124;03m\"\"\"Split a pathname.  Returns tuple \"(head, tail)\" where \"tail\" is\u001b[39;00m\n\u001b[1;32m    102\u001b[0m \u001b[38;5;124;03m    everything after the final slash.  Either part may be empty.\"\"\"\u001b[39;00m\n\u001b[0;32m--> 103\u001b[0m     p \u001b[38;5;241m=\u001b[39m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfspath\u001b[49m\u001b[43m(\u001b[49m\u001b[43mp\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    104\u001b[0m     sep \u001b[38;5;241m=\u001b[39m _get_sep(p)\n\u001b[1;32m    105\u001b[0m     i \u001b[38;5;241m=\u001b[39m p\u001b[38;5;241m.\u001b[39mrfind(sep) \u001b[38;5;241m+\u001b[39m \u001b[38;5;241m1\u001b[39m\n",
      "\u001b[0;31mTypeError\u001b[0m: expected str, bytes or os.PathLike object, not NoneType"
     ]
    }
   ],
   "source": [
    "opt = init()\n",
    "opt.figpath = figpath\n",
    "\n",
    "c = Caustique(opt)\n",
    "# compute just once\n",
    "z = c.wave()\n",
    "\n",
    "for H_ in c.opt.H*np.logspace(-1, 1, N_scan, base=base):\n",
    "    opt = init()\n",
    "    opt.figpath = figpath\n",
    "    c = Caustique(opt)\n",
    "\n",
    "    print(f'H = {H_:.3f}')\n",
    "    c.opt.H = H_\n",
    "    c.opt.tag = f'{opt.tag}_H_{H_:.3f}'\n",
    "    videoname = f'{opt.figpath}/{c.opt.tag}.{opt.vext}'\n",
    "    if not os.path.isfile(videoname):\n",
    "        url=c.plot(z, videoname=videoname)\n",
    "    c.show(videoname)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:59.431535Z",
     "iopub.status.busy": "2022-07-28T16:28:59.430866Z",
     "iopub.status.idle": "2022-07-28T16:28:59.474115Z",
     "shell.execute_reply": "2022-07-28T16:28:59.474566Z"
    }
   },
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "expected str, bytes or os.PathLike object, not NoneType",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Input \u001b[0;32mIn [14]\u001b[0m, in \u001b[0;36m<cell line: 4>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m opt \u001b[38;5;241m=\u001b[39m init()\n\u001b[1;32m      2\u001b[0m opt\u001b[38;5;241m.\u001b[39mfigpath \u001b[38;5;241m=\u001b[39m figpath\n\u001b[0;32m----> 4\u001b[0m c \u001b[38;5;241m=\u001b[39m \u001b[43mCaustique\u001b[49m\u001b[43m(\u001b[49m\u001b[43mopt\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m      5\u001b[0m z \u001b[38;5;241m=\u001b[39m c\u001b[38;5;241m.\u001b[39mwave()\n\u001b[1;32m      7\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m variation_ \u001b[38;5;129;01min\u001b[39;00m np\u001b[38;5;241m.\u001b[39mlogspace(\u001b[38;5;241m-\u001b[39m\u001b[38;5;241m2\u001b[39m, \u001b[38;5;241m0\u001b[39m, N_scan, base\u001b[38;5;241m=\u001b[39m\u001b[38;5;241m10\u001b[39m, endpoint\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m):\n",
      "Input \u001b[0;32mIn [4]\u001b[0m, in \u001b[0;36mCaustique.__init__\u001b[0;34m(self, opt)\u001b[0m\n\u001b[1;32m     86\u001b[0m \u001b[38;5;66;03m# https://stackoverflow.com/questions/16878315/what-is-the-right-way-to-treat-python-argparse-namespace-as-a-dictionary\u001b[39;00m\n\u001b[1;32m     87\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39md \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mvars\u001b[39m(opt)\n\u001b[0;32m---> 88\u001b[0m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mmakedirs\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mopt\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfigpath\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mexist_ok\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;28;43;01mTrue\u001b[39;49;00m\u001b[43m)\u001b[49m\n\u001b[1;32m     89\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath \u001b[38;5;241m=\u001b[39m os\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39mjoin(\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m/tmp\u001b[39m\u001b[38;5;124m'\u001b[39m, \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mopt\u001b[38;5;241m.\u001b[39mfigpath)\n\u001b[1;32m     90\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m opt\u001b[38;5;241m.\u001b[39mverbose: \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00m\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath\u001b[38;5;132;01m=}\u001b[39;00m\u001b[38;5;124m'\u001b[39m)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/os.py:210\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n\u001b[1;32m    200\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mmakedirs\u001b[39m(name, mode\u001b[38;5;241m=\u001b[39m\u001b[38;5;241m0o777\u001b[39m, exist_ok\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m):\n\u001b[1;32m    201\u001b[0m     \u001b[38;5;124;03m\"\"\"makedirs(name [, mode=0o777][, exist_ok=False])\u001b[39;00m\n\u001b[1;32m    202\u001b[0m \n\u001b[1;32m    203\u001b[0m \u001b[38;5;124;03m    Super-mkdir; create a leaf directory and all intermediate ones.  Works like\u001b[39;00m\n\u001b[0;32m   (...)\u001b[0m\n\u001b[1;32m    208\u001b[0m \n\u001b[1;32m    209\u001b[0m \u001b[38;5;124;03m    \"\"\"\u001b[39;00m\n\u001b[0;32m--> 210\u001b[0m     head, tail \u001b[38;5;241m=\u001b[39m \u001b[43mpath\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msplit\u001b[49m\u001b[43m(\u001b[49m\u001b[43mname\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    211\u001b[0m     \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m tail:\n\u001b[1;32m    212\u001b[0m         head, tail \u001b[38;5;241m=\u001b[39m path\u001b[38;5;241m.\u001b[39msplit(head)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/posixpath.py:103\u001b[0m, in \u001b[0;36msplit\u001b[0;34m(p)\u001b[0m\n\u001b[1;32m    100\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21msplit\u001b[39m(p):\n\u001b[1;32m    101\u001b[0m     \u001b[38;5;124;03m\"\"\"Split a pathname.  Returns tuple \"(head, tail)\" where \"tail\" is\u001b[39;00m\n\u001b[1;32m    102\u001b[0m \u001b[38;5;124;03m    everything after the final slash.  Either part may be empty.\"\"\"\u001b[39;00m\n\u001b[0;32m--> 103\u001b[0m     p \u001b[38;5;241m=\u001b[39m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfspath\u001b[49m\u001b[43m(\u001b[49m\u001b[43mp\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    104\u001b[0m     sep \u001b[38;5;241m=\u001b[39m _get_sep(p)\n\u001b[1;32m    105\u001b[0m     i \u001b[38;5;241m=\u001b[39m p\u001b[38;5;241m.\u001b[39mrfind(sep) \u001b[38;5;241m+\u001b[39m \u001b[38;5;241m1\u001b[39m\n",
      "\u001b[0;31mTypeError\u001b[0m: expected str, bytes or os.PathLike object, not NoneType"
     ]
    }
   ],
   "source": [
    "opt = init()\n",
    "opt.figpath = figpath\n",
    "\n",
    "c = Caustique(opt)\n",
    "z = c.wave()\n",
    "\n",
    "for variation_ in np.logspace(-2, 0, N_scan, base=10, endpoint=False):\n",
    "    opt = init()\n",
    "    opt.figpath = figpath\n",
    "    c = Caustique(opt)\n",
    "    print(f'variation = {variation_:.3f}')\n",
    "    c.opt.variation = variation_\n",
    "    c.opt.tag = f'{opt.tag}_variation_{variation_:.3f}'\n",
    "    videoname = f'{opt.figpath}/{c.opt.tag}.{opt.vext}'\n",
    "    if not os.path.isfile(videoname):\n",
    "        url=c.plot(z, videoname=videoname)\n",
    "    c.show(videoname)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-28T16:28:59.481143Z",
     "iopub.status.busy": "2022-07-28T16:28:59.480560Z",
     "iopub.status.idle": "2022-07-28T16:28:59.528798Z",
     "shell.execute_reply": "2022-07-28T16:28:59.529270Z"
    },
    "tags": []
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "======gamma======\n"
     ]
    },
    {
     "ename": "TypeError",
     "evalue": "expected str, bytes or os.PathLike object, not NoneType",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Input \u001b[0;32mIn [15]\u001b[0m, in \u001b[0;36m<cell line: 1>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      4\u001b[0m opt \u001b[38;5;241m=\u001b[39m init()\n\u001b[1;32m      5\u001b[0m opt\u001b[38;5;241m.\u001b[39mfigpath \u001b[38;5;241m=\u001b[39m figpath\n\u001b[0;32m----> 7\u001b[0m c \u001b[38;5;241m=\u001b[39m \u001b[43mCaustique\u001b[49m\u001b[43m(\u001b[49m\u001b[43mopt\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m      8\u001b[0m c\u001b[38;5;241m.\u001b[39md[variable] \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m=\u001b[39m modul\n\u001b[1;32m      9\u001b[0m c\u001b[38;5;241m.\u001b[39mopt\u001b[38;5;241m.\u001b[39mtag \u001b[38;5;241m=\u001b[39m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mopt\u001b[38;5;241m.\u001b[39mtag\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m_\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mvariable\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m_modul_\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mmodul\u001b[38;5;132;01m:\u001b[39;00m\u001b[38;5;124m.3f\u001b[39m\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m'\u001b[39m\n",
      "Input \u001b[0;32mIn [4]\u001b[0m, in \u001b[0;36mCaustique.__init__\u001b[0;34m(self, opt)\u001b[0m\n\u001b[1;32m     86\u001b[0m \u001b[38;5;66;03m# https://stackoverflow.com/questions/16878315/what-is-the-right-way-to-treat-python-argparse-namespace-as-a-dictionary\u001b[39;00m\n\u001b[1;32m     87\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39md \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mvars\u001b[39m(opt)\n\u001b[0;32m---> 88\u001b[0m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mmakedirs\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mopt\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfigpath\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mexist_ok\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;28;43;01mTrue\u001b[39;49;00m\u001b[43m)\u001b[49m\n\u001b[1;32m     89\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath \u001b[38;5;241m=\u001b[39m os\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39mjoin(\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m/tmp\u001b[39m\u001b[38;5;124m'\u001b[39m, \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mopt\u001b[38;5;241m.\u001b[39mfigpath)\n\u001b[1;32m     90\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m opt\u001b[38;5;241m.\u001b[39mverbose: \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124mf\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;132;01m{\u001b[39;00m\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcachepath\u001b[38;5;132;01m=}\u001b[39;00m\u001b[38;5;124m'\u001b[39m)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/os.py:210\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n\u001b[1;32m    200\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mmakedirs\u001b[39m(name, mode\u001b[38;5;241m=\u001b[39m\u001b[38;5;241m0o777\u001b[39m, exist_ok\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m):\n\u001b[1;32m    201\u001b[0m     \u001b[38;5;124;03m\"\"\"makedirs(name [, mode=0o777][, exist_ok=False])\u001b[39;00m\n\u001b[1;32m    202\u001b[0m \n\u001b[1;32m    203\u001b[0m \u001b[38;5;124;03m    Super-mkdir; create a leaf directory and all intermediate ones.  Works like\u001b[39;00m\n\u001b[0;32m   (...)\u001b[0m\n\u001b[1;32m    208\u001b[0m \n\u001b[1;32m    209\u001b[0m \u001b[38;5;124;03m    \"\"\"\u001b[39;00m\n\u001b[0;32m--> 210\u001b[0m     head, tail \u001b[38;5;241m=\u001b[39m \u001b[43mpath\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msplit\u001b[49m\u001b[43m(\u001b[49m\u001b[43mname\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    211\u001b[0m     \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m tail:\n\u001b[1;32m    212\u001b[0m         head, tail \u001b[38;5;241m=\u001b[39m path\u001b[38;5;241m.\u001b[39msplit(head)\n",
      "File \u001b[0;32m/usr/local/Cellar/python@3.9/3.9.13_1/Frameworks/Python.framework/Versions/3.9/lib/python3.9/posixpath.py:103\u001b[0m, in \u001b[0;36msplit\u001b[0;34m(p)\u001b[0m\n\u001b[1;32m    100\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21msplit\u001b[39m(p):\n\u001b[1;32m    101\u001b[0m     \u001b[38;5;124;03m\"\"\"Split a pathname.  Returns tuple \"(head, tail)\" where \"tail\" is\u001b[39;00m\n\u001b[1;32m    102\u001b[0m \u001b[38;5;124;03m    everything after the final slash.  Either part may be empty.\"\"\"\u001b[39;00m\n\u001b[0;32m--> 103\u001b[0m     p \u001b[38;5;241m=\u001b[39m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfspath\u001b[49m\u001b[43m(\u001b[49m\u001b[43mp\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    104\u001b[0m     sep \u001b[38;5;241m=\u001b[39m _get_sep(p)\n\u001b[1;32m    105\u001b[0m     i \u001b[38;5;241m=\u001b[39m p\u001b[38;5;241m.\u001b[39mrfind(sep) \u001b[38;5;241m+\u001b[39m \u001b[38;5;241m1\u001b[39m\n",
      "\u001b[0;31mTypeError\u001b[0m: expected str, bytes or os.PathLike object, not NoneType"
     ]
    }
   ],
   "source": [
    "for variable in ['gamma', 'sf_0', 'B_sf', 'B_theta', 'V_X', 'B_V']: #  'theta', 'V_Y'\n",
    "    print(f'======{variable}======')\n",
    "    for modul in np.logspace(-1, 1, N_scan, base=base):\n",
    "        opt = init()\n",
    "        opt.figpath = figpath\n",
    "\n",
    "        c = Caustique(opt)\n",
    "        c.d[variable] *= modul\n",
    "        c.opt.tag = f'{opt.tag}_{variable}_modul_{modul:.3f}'\n",
    "        videoname = f'{opt.figpath}/{c.opt.tag}.{opt.vext}'\n",
    "\n",
    "        print(f'{variable}={variable}(default)*{modul:.3f}={c.d[variable]:.3E}')\n",
    "        if not os.path.isfile(videoname):\n",
    "            print('Doing ', videoname)\n",
    "            z = c.wave()\n",
    "            mcname = f'{opt.figpath}/{c.opt.tag}'\n",
    "            if False: #not os.path.isfile(f'{mcname}{c.mc.vext}'): \n",
    "                print('Doing ', f'{mcname}{c.mc.vext}')\n",
    "                c.mc.anim_save(z.swapaxes(0, 1), f'{mcname}')\n",
    "            url=c.plot(z, videoname=videoname)\n",
    "        c.show(videoname)"
   ]
  }
 ],
 "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.9.13"
  },
  "nteract": {
   "version": "0.28.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}