Skip to content

Instantly share code, notes, and snippets.

@davidwhogg
Forked from adrn/ImagePhaseSandbox.ipynb
Created October 26, 2019 16:44
Show Gist options
  • Select an option

  • Save davidwhogg/4dcc600188772201df6fce2489a86a44 to your computer and use it in GitHub Desktop.

Select an option

Save davidwhogg/4dcc600188772201df6fce2489a86a44 to your computer and use it in GitHub Desktop.
ImagePhaseSandbox.ipynb
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "import astropy.coordinates as coord\nimport astropy.units as u\nimport matplotlib as mpl\nimport matplotlib.pyplot as plt\n%matplotlib inline\nimport numpy as np\n\n# from scipy.fftpack import rfft as fft, irfft as ifft\nfrom scipy.fftpack import fft2 as fft, ifft2 as ifft\nfrom PIL import Image as pil",
"execution_count": 1,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "# import cv2\n# vidcap = cv2.VideoCapture('/Users/apricewhelan/Downloads/48728220883_b87ea8cf30_vm.mp4')\n# success, image = vidcap.read()\n# count = 0\n# while success:\n# if count == 0 or count == 31:\n# cv2.imwrite(\"frame%d.jpg\" % count, image) # save frame as JPEG file \n# success, image = vidcap.read()\n# count += 1",
"execution_count": 2,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "rgbs = []\nfor filename in ['frame0.jpg', 'frame30.jpg']:\n im = pil.open(filename)\n rgb = np.array(im).astype(np.float64)\n rgbs.append(rgb)",
"execution_count": 3,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "# for rgb in rgbs:\n# fig, ax = plt.subplots(figsize=(10, 10))\n# ax.imshow(rgb[..., 0], cmap='Greys')\n# ax.set_aspect('equal')",
"execution_count": 35,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "all_ffts = []\nfor rgb in rgbs:\n for band in range(3):\n all_ffts.append(fft(rgb[..., band]))\n \n# fig, ax = plt.subplots(figsize=(6, 6))\n# ax.imshow(ifft(fft(rgb)), cmap='Greys')\n# ax.set_aspect('equal')",
"execution_count": 5,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "def get_ims(data1, data2, nslice, save_dont_return=False):\n f1 = np.stack((data1.real, data1.imag), axis=0)\n f2 = np.stack((data2.real, data2.imag), axis=0)\n \n amp1 = np.abs(data1)\n amp2 = np.abs(data2)\n \n uf1 = f1 / amp1[None]\n uf2 = f2 / amp2[None]\n\n cross = np.cross(uf1, uf2, axis=0)\n theta2 = np.arcsin(cross)\n thetas = np.linspace(0, theta2, nslice)\n dtheta = thetas[1] - thetas[0]\n \n frames = np.concatenate((np.arange(-4*nslice, 0), \n np.arange(0, nslice),\n np.arange(nslice, nslice + 4*nslice)))\n \n ims = [] \n for j, i in enumerate(frames):\n theta = dtheta * i\n if i < 0:\n fac = amp1\n elif i > nslice:\n fac = amp2\n else:\n fac = (amp2 - amp1) / nslice * i + amp1\n \n R = np.array([[np.cos(theta), np.sin(theta)],\n [-np.sin(theta), np.cos(theta)]])\n rot_f1 = np.einsum('ijnm,inm->jnm', R, uf1)\n \n rot_f1 *= fac\n rot_f1 = rot_f1[0] + 1j*rot_f1[1]\n rot_f1 = ifft(rot_f1).real\n \n if save_dont_return:\n fig, ax = plt.subplots(figsize=(10, 8))\n ax.imshow(rot_f1.real, cmap='Greys')\n ax.set_aspect('equal')\n ax.xaxis.set_visible(False)\n ax.yaxis.set_visible(False)\n fig.tight_layout()\n fig.savefig(f'/Users/apricewhelan/projects/scratch/supernova-phase-bump/sn_frame_{j:03d}.png', \n dpi=150)\n plt.close(fig)\n \n else:\n ims.append(rot_f1)\n \n if not save_dont_return:\n return np.array(ims)",
"execution_count": 66,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "import glob, os\nfor filename in glob.glob('/Users/apricewhelan/projects/scratch/supernova-phase-bump/*'):\n num = int(filename.split('.')[0].split('_')[-1])\n new_filename = '/Users/apricewhelan/projects/scratch/supernova-phase-bump/sn_frame_{:03d}.png'.format(num)\n os.system('mv {} {}'.format(filename, new_filename))",
"execution_count": 73,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "nslice = 16\nims = get_ims(all_ffts[0], all_ffts[3], nslice, save_dont_return=True)",
"execution_count": 67,
"outputs": []
},
{
"metadata": {
"scrolled": false,
"trusted": true
},
"cell_type": "code",
"source": "# for im in [ims[0], ims[nslice], ims[nslice+nslice], ims[-1]]:\n# fig, ax = plt.subplots(figsize=(10, 8))\n# ax.imshow(im.real, cmap='Greys')\n# ax.set_aspect('equal')\n# ax.xaxis.set_visible(False)\n# ax.yaxis.set_visible(False)\n# fig.tight_layout()",
"execution_count": null,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "",
"execution_count": null,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "",
"execution_count": null,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "",
"execution_count": null,
"outputs": []
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "# max_theta = np.arcsin(cross[2])\n\nR = np.array([[np.cos(theta), -np.sin(theta)],\n [np.sin(theta), np.cos(theta)]])\nR.shape",
"execution_count": 113,
"outputs": [
{
"data": {
"text/plain": "(2, 2, 1584, 1860)"
},
"execution_count": 113,
"metadata": {},
"output_type": "execute_result"
}
]
},
{
"metadata": {
"trusted": true
},
"cell_type": "code",
"source": "plt.figure(figsize=(10, 10))\nplt.imshow(ifft(amp1 * np.exp(1j * phase2)).real - rgbs[0][..., 0], cmap='Greys')\n# plt.colorbar()",
"execution_count": 90,
"outputs": [
{
"data": {
"text/plain": "<matplotlib.image.AxesImage at 0x1c2f09e9b0>"
},
"execution_count": 90,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment