JonesCalculus_15b.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Defining Classes for Jones vectors, matrices, and density matrices\n",
    "#### ProfHuster 20151123\n",
    "#### Based on http://spie.org/x32380.xml\n",
    "----"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from IPython.display import display, Math, Latex\n",
    "from sympy import *\n",
    "from sympy.utilities.misc import filldedent\n",
    "from sympy import init_printing\n",
    "def my_display(x):\n",
    "    display(Math(latex(x)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Define a classes for Jones Vectors, Matrices, and Density Matrices"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "class JVec(Matrix):\n",
    "    \"\"\"\n",
    "    Jones vectors for describing light with pure polarization\n",
    "\n",
    "    Jones vectors are 2 element complex column vectors describing the polarization\n",
    "    of light. The elements are [[Ex],[Ey]] where Ex and Ey are the complex\n",
    "    amplitudes of the electric field vector. The convention is that Ex is horizontal,\n",
    "    and Ey is vertical. \n",
    "    They are oftern normalized such that the incident light beam has\n",
    "    conj(Ex)*Ex + conj(Ey)Ey=1.\n",
    "\n",
    "    Attributes:\n",
    "        All Matrix attributes\n",
    "    \"\"\"\n",
    "    def __new__(cls, *args):\n",
    "        newobj = Matrix.__new__(cls, *args)\n",
    "        if newobj.shape != (2, 1):\n",
    "            raise TypeError(\"JVec: shape must be (2,1)\")\n",
    "        return newobj\n",
    "\n",
    "# define functions to return common polarizations\n",
    "def LHP():\n",
    "    return JVec([1,0])\n",
    "def LVP():\n",
    "    return JVec([0,1])\n",
    "def LP(theta):\n",
    "    return JVec([cos(theta),sin(theta)])\n",
    "def RCP():\n",
    "    return JVec([1/np.sqrt(2),1j/np.sqrt(2)])\n",
    "def LCP():\n",
    "    return JVec([1/np.sqrt(2),-1j/np.sqrt(2)])\n",
    "\n",
    "class JMat(Matrix):\n",
    "    \"\"\"\n",
    "    Jones matrices describing optical elements that change light's E fields\n",
    "    \n",
    "    Jones matrices are 2x2 complex matrices that operate on a Jones vector,\n",
    "    a complex, two element column vector.\n",
    "\n",
    "    Attributes:\n",
    "        All Matrix attributes\n",
    "    \"\"\"\n",
    "    def __new__(cls, *args):\n",
    "        newobj = Matrix.__new__(cls, *args)\n",
    "        if newobj.shape != (2, 2):\n",
    "            raise TypeError(\"JMat: shape must be (2,2)\")\n",
    "        return newobj\n",
    "\n",
    "def JM_LHP():\n",
    "    return JMat([[1,0],[0,0]])\n",
    "def JM_LVP():\n",
    "    return JMat([[0,0],[0,1]])\n",
    "def JM_QWP():\n",
    "    return JMat([[1,0],[0,-I]])\n",
    "def JM_HWP():\n",
    "    return JMat([[1,0],[0,-1]])\n",
    "def JM_WP(phi):\n",
    "    return JMat([[exp(I*phi/2),0],[0,exp(-I*phi/2)]])\n",
    "def JM_ROT(theta):\n",
    "    return JMat([[cos(theta),sin(theta)],[-sin(theta),cos(theta)]])\n",
    "\n",
    "class JDensity(Matrix):\n",
    "    \"\"\"\n",
    "    A Jones density matrix is an extension of the Jones vectors to include partially\n",
    "    polarized light.\n",
    "    \n",
    "    Jones density matrices are 2x2 complex matrices\n",
    "\n",
    "    Attributes:\n",
    "        All Matrix attributes\n",
    "    \"\"\"\n",
    "    def __new__(cls, *args):\n",
    "        newobj = Matrix.__new__(cls, *args)\n",
    "        if newobj.shape != (2, 2):\n",
    "            raise TypeError(\"JDensity: shape must be (2,2)\")\n",
    "        return newobj"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "LHP, LVP, RCP, LCP, LP(alpha) =\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}1\\\\0\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}0\\\\1\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}0.707106781186547\\\\0.707106781186547 i\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}0.707106781186547\\\\- 0.707106781186547 i\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}\\cos{\\left (\\alpha \\right )}\\\\\\sin{\\left (\\alpha \\right )}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "JM_LHP, JM_LVP, JM_QWP, JM_HWP, JM_ROT(beta) = \n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}1 & 0\\\\0 & 0\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}0 & 0\\\\0 & 1\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}1 & 0\\\\0 & - i\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}1 & 0\\\\0 & -1\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}\\cos{\\left (\\beta \\right )} & \\sin{\\left (\\beta \\right )}\\\\- \\sin{\\left (\\beta \\right )} & \\cos{\\left (\\beta \\right )}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Wave Plate =\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}e^{\\frac{i \\phi}{2}} & 0\\\\0 & e^{- \\frac{i \\phi}{2}}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Rotated WP = \n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}e^{\\frac{i \\phi}{2}} \\cos^{2}{\\left (\\theta \\right )} + e^{- \\frac{i \\phi}{2}} \\sin^{2}{\\left (\\theta \\right )} & e^{\\frac{i \\phi}{2}} \\sin{\\left (\\theta \\right )} \\cos{\\left (\\theta \\right )} - e^{- \\frac{i \\phi}{2}} \\sin{\\left (\\theta \\right )} \\cos{\\left (\\theta \\right )}\\\\e^{\\frac{i \\phi}{2}} \\sin{\\left (\\theta \\right )} \\cos{\\left (\\theta \\right )} - e^{- \\frac{i \\phi}{2}} \\sin{\\left (\\theta \\right )} \\cos{\\left (\\theta \\right )} & e^{\\frac{i \\phi}{2}} \\sin^{2}{\\left (\\theta \\right )} + e^{- \\frac{i \\phi}{2}} \\cos^{2}{\\left (\\theta \\right )}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "0,0=\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "$$e^{\\frac{i \\phi}{2}} \\cos^{2}{\\left (\\theta \\right )} + e^{- \\frac{i \\phi}{2}} \\sin^{2}{\\left (\\theta \\right )}$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha = Symbol(\"alpha\", real=True)\n",
    "A = LHP()\n",
    "B = LVP()\n",
    "C = RCP()\n",
    "D = LCP()\n",
    "E = LP(alpha)\n",
    "print \"LHP, LVP, RCP, LCP, LP(alpha) =\"\n",
    "my_display(A)\n",
    "my_display(B)\n",
    "my_display(C)\n",
    "my_display(D)\n",
    "my_display(E)\n",
    "\n",
    "beta = Symbol('beta', real=True)\n",
    "JMA = JM_LHP()\n",
    "JMB = JM_LVP()\n",
    "JMQW = JM_QWP()\n",
    "JMHW = JM_HWP()\n",
    "JMRotBeta = JM_ROT(beta)\n",
    "print \"JM_LHP, JM_LVP, JM_QWP, JM_HWP, JM_ROT(beta) = \"\n",
    "my_display(JMA)\n",
    "my_display(JMB)\n",
    "my_display(JMQW)\n",
    "my_display(JMHW)\n",
    "my_display(JMRotBeta)\n",
    "\n",
    "# arb wave plate\n",
    "print \"\\nWave Plate =\"\n",
    "phi = Symbol(\"phi\", real=True)\n",
    "my_display(JM_WP(phi))\n",
    "# rotated arb wave plate\n",
    "print \"\\nRotated WP = \"\n",
    "theta = Symbol(\"theta\", real=True)\n",
    "JM_E = JM_ROT(-theta) * JM_WP(phi) * JM_ROT(theta)\n",
    "my_display(JM_E)\n",
    "print \"\\n0,0=\"\n",
    "my_display(JM_E[0,0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}1\\\\0\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "my_display(LHP())\n",
    "\n",
    "# arb wave plate\n",
    "#my_display(JM_WP(phi))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "0,0=\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "$$\\left(- i \\sin{\\left (\\frac{\\phi}{2} \\right )} + \\cos{\\left (\\frac{\\phi}{2} \\right )}\\right) \\sin^{2}{\\left (\\theta \\right )} + \\left(i \\sin{\\left (\\frac{\\phi}{2} \\right )} + \\cos{\\left (\\frac{\\phi}{2} \\right )}\\right) \\cos^{2}{\\left (\\theta \\right )}$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "====\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}\\left(- i \\sin{\\left (\\frac{\\phi}{2} \\right )} + \\cos{\\left (\\frac{\\phi}{2} \\right )}\\right) \\left(\\left(i \\sin{\\left (\\phi \\right )} + \\cos{\\left (\\phi \\right )}\\right) \\cos^{2}{\\left (\\theta \\right )} + \\sin^{2}{\\left (\\theta \\right )}\\right) & i \\sin{\\left (\\frac{\\phi}{2} \\right )} \\sin{\\left (2 \\theta \\right )}\\\\i \\sin{\\left (\\frac{\\phi}{2} \\right )} \\sin{\\left (2 \\theta \\right )} & \\left(- i \\sin{\\left (\\frac{\\phi}{2} \\right )} + \\cos{\\left (\\frac{\\phi}{2} \\right )}\\right) \\left(\\left(i \\sin{\\left (\\phi \\right )} + \\cos{\\left (\\phi \\right )}\\right) \\sin^{2}{\\left (\\theta \\right )} + \\cos^{2}{\\left (\\theta \\right )}\\right)\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "a = JM_E[0,0]\n",
    "print \"\\n0,0=\"\n",
    "b = a.rewrite(cos)\n",
    "my_display(b)\n",
    "b = a.rewrite(exp)\n",
    "print \"====\\n\"\n",
    "my_display(simplify(JM_E).rewrite(cos))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Define some basic Jones Matrices\n",
    "* Linear & Circular Polarizers\n",
    "* Quarter Wave Plates\n",
    "* Rotator"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "rcp=\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}\\frac{\\sqrt{2}}{2}\\\\\\frac{\\sqrt{2} i}{2}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}0\\\\1\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}0\\\\e^{\\frac{i \\phi}{2}} \\sin^{2}{\\left (\\beta \\right )} + e^{- \\frac{i \\phi}{2}} \\cos^{2}{\\left (\\beta \\right )}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\left(e^{\\frac{i \\phi}{2}} \\sin^{2}{\\left (\\beta \\right )} + e^{- \\frac{i \\phi}{2}} \\cos^{2}{\\left (\\beta \\right )}\\right) \\left(e^{\\frac{i \\overline{\\phi}}{2}} \\cos^{2}{\\left (\\overline{\\beta} \\right )} + e^{- \\frac{i \\overline{\\phi}}{2}} \\sin^{2}{\\left (\\overline{\\beta} \\right )}\\right)$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Polarization states\n",
    "lvp = Matrix([[0], [1]])\n",
    "rcp = (1/sqrt(2)) * Matrix([[1],[I]])\n",
    "print \"rcp=\"\n",
    "my_display(rcp)\n",
    "# Devices\n",
    "LHP = Matrix([[1, 0],[0, 0]])\n",
    "LVP = Matrix([[0,0],[0,1]])\n",
    "# Next a linear polarizer at an angle theta from the x-axis\n",
    "theta = Symbol('theta')\n",
    "LP_Theta = Matrix([[cos(theta)**2, sin(theta)*cos(theta)], \\\n",
    "                   [sin(theta)*cos(theta), sin(theta)**2]])\n",
    "\n",
    "# right and left circular polarizers\n",
    "half = Rational(1,2)\n",
    "RCP = Matrix([[half, I/2],[-I/2,half]])\n",
    "LCP = Matrix([[half, -I/2],[I/2,half]])\n",
    "\n",
    "# retarders\n",
    "phi = Symbol('phi')\n",
    "WP_phi = Matrix([[exp(I*phi/2), 0], [0, exp(-I*phi/2)]])\n",
    "QW_SH = Matrix([[1, 0], [0, -I]])\n",
    "QW_SV = Matrix([[1, 0], [0,  I]])\n",
    "beta = Symbol('beta')\n",
    "R_Beta = Matrix([[cos(beta), -sin(beta)],[sin(beta), cos(beta)]])\n",
    "R_IBeta = Matrix([[cos(beta), sin(beta)],[-sin(beta), cos(beta)]])\n",
    "my_display(lvp)\n",
    "exitState = LVP * R_IBeta * WP_phi * R_Beta * lvp\n",
    "my_display(exitState)\n",
    "exitI = exitState.H * exitState\n",
    "my_display(exitI[0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Now plot versus angle"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left(e^{- \\frac{i \\pi}{8}} \\sin^{2}{\\left (\\overline{\\beta} \\right )} + e^{\\frac{i \\pi}{8}} \\cos^{2}{\\left (\\overline{\\beta} \\right )}\\right) \\left(e^{\\frac{i \\pi}{8}} \\sin^{2}{\\left (\\beta \\right )} + e^{- \\frac{i \\pi}{8}} \\cos^{2}{\\left (\\beta \\right )}\\right)$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sympy.plotting import plot\n",
    "I_Plot = exitI[0].subs([(phi,pi/4)])\n",
    "plot(I_Plot,(beta,0,6.28))\n",
    "my_display(I_Plot)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Define some Jones Density Matrices for partially polarized light"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}\\frac{1}{2} & - 0.4 i\\\\0.4 i & \\frac{1}{2}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Non-polarized light\n",
    "rhoNP = Matrix([[half, 0], [0, half]])\n",
    "# some pure polarizations\n",
    "rhoLHP = Matrix([[1,0],[0,0]])\n",
    "rhoLVP = Matrix([[0,0],[0,1]])\n",
    "rhoRCP = Matrix([[half,I/2],[-I/2,half]])\n",
    "rhoLCP = Matrix([[half,-I/2],[I/2,half]])\n",
    "# partial Linear and Circular polarizations\n",
    "p = Symbol('p')\n",
    "rhoPartLP = Matrix([[p,0],[0,(1-p)]])\n",
    "rhoPartCP = Matrix([[half,I*(p-half)],[-I*(p-half),half]]) \n",
    "my_display(rhoPartCP.subs(p,.1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Define the system for measuring polarization state\n",
    "The system is a rotating quarter wave plate followed by a linear polarizer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}- i \\sin^{2}{\\left (\\beta \\right )} + \\cos^{2}{\\left (\\beta \\right )} & \\frac{1}{2} \\left(1 + i\\right) \\sin{\\left (2 \\beta \\right )}\\\\0 & 0\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "M_Sys = simplify(factor(LHP * R_Beta * QW_SH * R_IBeta))\n",
    "my_display(M_Sys)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}\\left(\\frac{1}{2} + \\frac{i}{2}\\right) \\cos{\\left (2 \\beta \\right )} + \\frac{1}{2} - \\frac{i}{2} & \\left(\\frac{1}{2} + \\frac{i}{2}\\right) \\sin{\\left (2 \\beta \\right )}\\\\0 & 0\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "M_Sys1 = Matrix([(half*(1-I)+half*(1+I)*cos(2*beta),(1+I)*half*sin(2*beta)),(0,0)])\n",
    "my_display(M_Sys1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Now define a general density matrix\n",
    "The density matrix has to be Hermitian, so for a normalized beam (I = 1), there are only three parameters:\n",
    "\n",
    "\n",
    "$\\left( \\begin{array}{ccc}\n",
    "a & b + i c \\\\\n",
    "b - i c & (1-a) \\end{array} \\right)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}a & b + i c\\\\b - i c & - a + 1\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$4 a^{2} - 4 a + 4 b^{2} + 4 c^{2} + 1$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$$\\frac{i a}{2} \\cos{\\left (2 \\beta \\right )} + \\frac{a}{2} \\cos{\\left (2 \\left(\\beta + \\overline{\\beta}\\right) \\right )} - \\frac{i a}{2} \\cos{\\left (2 \\overline{\\beta} \\right )} + \\frac{a}{2} + \\frac{i b}{2} \\sin{\\left (2 \\beta \\right )} + \\frac{b}{2} \\sin{\\left (2 \\left(\\beta + \\overline{\\beta}\\right) \\right )} - \\frac{i b}{2} \\sin{\\left (2 \\overline{\\beta} \\right )} + \\frac{c}{2} \\sin{\\left (2 \\beta \\right )} - \\frac{i c}{2} \\sin{\\left (2 \\left(\\beta - \\overline{\\beta}\\right) \\right )} + \\frac{c}{2} \\sin{\\left (2 \\overline{\\beta} \\right )} + \\frac{1}{4} \\cos{\\left (2 \\left(\\beta - \\overline{\\beta}\\right) \\right )} - \\frac{1}{4} \\cos{\\left (2 \\left(\\beta + \\overline{\\beta}\\right) \\right )}$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "a = Symbol('a', real=True)\n",
    "b = Symbol('b', real=True)\n",
    "c = Symbol('c', real=True)\n",
    "rhoGen = Matrix([[a, b+I*c], [b-I*c, (1-a)]])\n",
    "my_display(rhoGen)\n",
    "\n",
    "Pol = 2 * trace(rhoGen**2) / (trace(rhoGen)**2) - 1\n",
    "my_display(simplify(Pol))\n",
    "#\n",
    "M_Simp = simplify(factor(M_Sys))\n",
    "M_2 = M_Simp.H * M_Simp\n",
    "I_Rot = simplify(factor(trace(simplify(factor(rhoGen * M_2)))))\n",
    "my_display(I_Rot)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$0.05 \\sin{\\left (2 \\beta \\right )} + 0.05 i \\sin{\\left (2 \\beta \\right )} - 0.05 i \\sin{\\left (2 \\left(\\beta - \\overline{\\beta}\\right) \\right )} + 0.05 \\sin{\\left (2 \\left(\\beta + \\overline{\\beta}\\right) \\right )} + 0.05 \\sin{\\left (2 \\overline{\\beta} \\right )} - 0.05 i \\sin{\\left (2 \\overline{\\beta} \\right )} + 0.25 i \\cos{\\left (2 \\beta \\right )} + \\frac{1}{4} \\cos{\\left (2 \\left(\\beta - \\overline{\\beta}\\right) \\right )} - 0.25 i \\cos{\\left (2 \\overline{\\beta} \\right )} + 0.25$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sympy.plotting import plot\n",
    "I_Plot = I_Rot.subs([(a,.5),(b,0.1),(c,0.1)])\n",
    "plot(I_Plot,(beta,0,6.28))\n",
    "my_display(I_Plot)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}a & b + i c\\\\b - i c & - a + 1\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "my_display(rhoGen)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}