beam2stl.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import stl\n",
    "from stl import mesh\n",
    "import numpy as np\n",
    "import mwa_hyperbeam\n",
    "\n",
    "# beam = mwa_hyperbeam.FEEBeam(\"/Users/dev/Library/CloudStorage/OneDrive-Curtin/Code/mwa_hyperbeam/mwa_full_embedded_element_pattern.h5\")\n",
    "beam = mwa_hyperbeam.AnalyticBeam(bowties_per_row=8)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "min=-97.0765464692417, max=<built-in function max>\n"
     ]
    }
   ],
   "source": [
    "from stl import mesh\n",
    "\n",
    "radius = 1.0\n",
    "resolution = 512\n",
    "\n",
    "# Create a sphere mesh\n",
    "az = np.linspace(0, 2 * np.pi, resolution)\n",
    "za = np.linspace(0, np.pi / 2, resolution)\n",
    "za, az = np.meshgrid(za, az)\n",
    "az, za = az.flatten(), za.flatten()\n",
    "\n",
    "freq = 154880000\n",
    "# ndips = 16\n",
    "ndips = 64\n",
    "\n",
    "delays = [0] * ndips\n",
    "amps = [1.0] * ndips\n",
    "norm_to_zenith = True\n",
    "# array_latitude_rad = None\n",
    "array_latitude_rad = -0.4660608448386394 # MWA latitude\n",
    "iau_order = False\n",
    "# jones = beam.calc_jones_array(\n",
    "#     az, za, freq, delays, amps, norm_to_zenith, array_latitude_rad, iau_order)\n",
    "jones = beam.calc_jones_array(\n",
    "    az, za, freq, delays, amps, array_latitude_rad, norm_to_zenith)\n",
    "\n",
    "\n",
    "mag_x = np.abs(jones[:, 0])**2 + np.abs(jones[:, 1]) ** 2\n",
    "mag_x = np.log(mag_x)\n",
    "min = np.min(mag_x)\n",
    "print(f\"{min=}, {max=}\")\n",
    "max = np.max(mag_x)\n",
    "norm_x = (mag_x - min) / (max - min)\n",
    "norm_x **= 2\n",
    "norm_x = np.nan_to_num(norm_x, nan=0)\n",
    "\n",
    "\n",
    "x = norm_x * np.sin(za) * np.cos(az)\n",
    "y = norm_x * np.sin(za) * np.sin(az)\n",
    "z = norm_x * np.cos(za)\n",
    "\n",
    "# Create the vertices and faces of the mesh\n",
    "vertices = np.column_stack((x.ravel(), y.ravel(), z.ravel()))\n",
    "faces = []\n",
    "for i in range(resolution - 1):\n",
    "    for j in range(resolution - 1):\n",
    "        v0 = i * resolution + j\n",
    "        v1 = v0 + 1\n",
    "        v2 = (i + 1) * resolution + j\n",
    "        v3 = v2 + 1\n",
    "        faces.extend([(v0, v1, v2), (v1, v3, v2)])\n",
    "\n",
    "# Create the mesh object\n",
    "sphere_mesh = mesh.Mesh(np.zeros(faces.__len__(), dtype=mesh.Mesh.dtype))\n",
    "for i, f in enumerate(faces):\n",
    "    for j in range(3):\n",
    "        sphere_mesh.vectors[i][j] = vertices[f[j], :]\n",
    "sphere_mesh.rotate([1.0, 0.0, 0.0], np.radians(90))\n",
    "\n",
    "# Save the mesh to an STL file\n",
    "# sphere_mesh.save(f'../public/objects/beam_f{freq:09d}.stl')\n",
    "sphere_mesh.save(f'beam_f{freq:09d}.stl')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Figure size 640x480 with 0 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# from mpl_toolkits import mplot3d\n",
    "# import matplotlib.pyplot as plt\n",
    "\n",
    "# # plot sphere mesh with a 3d projection\n",
    "# figure = plt.figure()\n",
    "# axes = mplot3d.Axes3D(figure, pro)\n",
    "# axes.add_collection3d(mplot3d.art3d.Poly3DCollection(sphere_mesh.vectors))\n",
    "# plt.show()\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}