|  | 
|  | 1 | +{ | 
|  | 2 | + "cells": [ | 
|  | 3 | +  { | 
|  | 4 | +   "cell_type": "code", | 
|  | 5 | +   "execution_count": null, | 
|  | 6 | +   "id": "e375f8cb", | 
|  | 7 | +   "metadata": {}, | 
|  | 8 | +   "outputs": [], | 
|  | 9 | +   "source": [ | 
|  | 10 | +    "import pygrc as gr\n", | 
|  | 11 | +    "import math\n", | 
|  | 12 | +    "import numpy as np\n", | 
|  | 13 | +    "import matplotlib.pyplot as plt" | 
|  | 14 | +   ] | 
|  | 15 | +  }, | 
|  | 16 | +  { | 
|  | 17 | +   "cell_type": "code", | 
|  | 18 | +   "execution_count": null, | 
|  | 19 | +   "id": "8332c16c", | 
|  | 20 | +   "metadata": {}, | 
|  | 21 | +   "outputs": [], | 
|  | 22 | +   "source": [ | 
|  | 23 | +    "galaxy = [\n", | 
|  | 24 | +    "\"NGC4217\",\n", | 
|  | 25 | +    "\"NGC5055\",\n", | 
|  | 26 | +    "\"NGC5585\",\n", | 
|  | 27 | +    "\"NGC6015\",\n", | 
|  | 28 | +    "\"NGC6503\",\n", | 
|  | 29 | +    "\"NGC7331\",\n", | 
|  | 30 | +    "]" | 
|  | 31 | +   ] | 
|  | 32 | +  }, | 
|  | 33 | +  { | 
|  | 34 | +   "cell_type": "code", | 
|  | 35 | +   "execution_count": null, | 
|  | 36 | +   "id": "a81f659b", | 
|  | 37 | +   "metadata": { | 
|  | 38 | +    "scrolled": false | 
|  | 39 | +   }, | 
|  | 40 | +   "outputs": [], | 
|  | 41 | +   "source": [ | 
|  | 42 | +    "df_all = []\n", | 
|  | 43 | +    "for g in galaxy:\n", | 
|  | 44 | +    "    df = gr.Reader.read(filepath=\"../notebooks/data/\"+g+\"_rotmod.dat\")\n", | 
|  | 45 | +    "    df_all.append(df)\n", | 
|  | 46 | +    "    gr.Plot().overlap(df,\"Rad\",[\"Vobs\",\"Vgas\",\"Vbul\",\"Vdisk\"],\"Distance (kpc)\",\"Velocity\",g)" | 
|  | 47 | +   ] | 
|  | 48 | +  }, | 
|  | 49 | +  { | 
|  | 50 | +   "cell_type": "code", | 
|  | 51 | +   "execution_count": null, | 
|  | 52 | +   "id": "149622af", | 
|  | 53 | +   "metadata": {}, | 
|  | 54 | +   "outputs": [], | 
|  | 55 | +   "source": [ | 
|  | 56 | +    "r = np.linspace(0.01,50,2000)\n", | 
|  | 57 | +    "def mass(r, M0, R0):\n", | 
|  | 58 | +    "#def mass(r, M0, rc, R0,beta):\n", | 
|  | 59 | +    "    #return M0*(np.sqrt(R0/rc)*(r/(r+rc)))**(3*beta)\n", | 
|  | 60 | +    "    return M0*(1- (1+(r/R0))*np.exp(-r/R0))\n", | 
|  | 61 | +    "\n", | 
|  | 62 | +    "def mond(r, M0, rc, R0,b, beta):\n", | 
|  | 63 | +    "    a = 1.2e-10\n", | 
|  | 64 | +    "    G = 4.300e-6 #parsecs \n", | 
|  | 65 | +    "    m = mass(r,M0, R0)\n", | 
|  | 66 | +    "    f = (G*m/r)*(1 + b*(1+(r/rc)))#*10e-3\n", | 
|  | 67 | +    "    return np.sqrt(f)*10e-3\n", | 
|  | 68 | +    "\n", | 
|  | 69 | +    "def newton(r, M0, rc, R0,beta):\n", | 
|  | 70 | +    "    a = 1.2e-10\n", | 
|  | 71 | +    "    G = 4.30e-6 #parsecs \n", | 
|  | 72 | +    "    m = mass(r,M0, R0)\n", | 
|  | 73 | +    "    f = G*m/r\n", | 
|  | 74 | +    "    #f = (G*m/r)*(1/np.sqrt(2)) * np.sqrt(1 + np.sqrt(1 + (r**4) * (2*a/(G*m))**2))\n", | 
|  | 75 | +    "    return np.sqrt(f)*10e-3" | 
|  | 76 | +   ] | 
|  | 77 | +  }, | 
|  | 78 | +  { | 
|  | 79 | +   "cell_type": "code", | 
|  | 80 | +   "execution_count": null, | 
|  | 81 | +   "id": "9ad18e72", | 
|  | 82 | +   "metadata": {}, | 
|  | 83 | +   "outputs": [], | 
|  | 84 | +   "source": [ | 
|  | 85 | +    "for df in df_all:\n", | 
|  | 86 | +    "    print(df['Rad'].max())\n", | 
|  | 87 | +    "    print(df['Rad'].min())" | 
|  | 88 | +   ] | 
|  | 89 | +  }, | 
|  | 90 | +  { | 
|  | 91 | +   "cell_type": "code", | 
|  | 92 | +   "execution_count": null, | 
|  | 93 | +   "id": "d49750d3", | 
|  | 94 | +   "metadata": {}, | 
|  | 95 | +   "outputs": [], | 
|  | 96 | +   "source": [ | 
|  | 97 | +    "m_mond=[]\n", | 
|  | 98 | +    "m_newton=[]\n", | 
|  | 99 | +    "\n", | 
|  | 100 | +    "for df in df_all:\n", | 
|  | 101 | +    "    m_mond.append(gr.Fit(df['Rad'],df['Vobs'],1.,1.,3,.35,1.))\n", | 
|  | 102 | +    "    m_newton.append(gr.Fit(df['Rad'],df['Vobs'],1.,1.,3,1.))" | 
|  | 103 | +   ] | 
|  | 104 | +  }, | 
|  | 105 | +  { | 
|  | 106 | +   "cell_type": "code", | 
|  | 107 | +   "execution_count": null, | 
|  | 108 | +   "id": "7e2450de", | 
|  | 109 | +   "metadata": {}, | 
|  | 110 | +   "outputs": [], | 
|  | 111 | +   "source": [ | 
|  | 112 | +    "fit_mond = [0]*len(m_mond)\n", | 
|  | 113 | +    "fit_newton = [0]*len(m_newton)\n", | 
|  | 114 | +    "for i in range(len(m_mond)):\n", | 
|  | 115 | +    "    fit_mond[i] = m_mond[i].fit_lsq(mond, [(1.,1e14),(1.,10.),(2.,5.),(0.1,2),(0.1,2)],df_all[i]['errV'],[False,False,True,True,True])\n", | 
|  | 116 | +    "    fit_newton[i] = m_newton[i].fit_lsq(newton, [(1.,1e14),(1.,10.),(1,10),(0.1,2)],df_all[i]['errV'],[False,True,True,True])" | 
|  | 117 | +   ] | 
|  | 118 | +  }, | 
|  | 119 | +  { | 
|  | 120 | +   "cell_type": "code", | 
|  | 121 | +   "execution_count": null, | 
|  | 122 | +   "id": "f770703a", | 
|  | 123 | +   "metadata": {}, | 
|  | 124 | +   "outputs": [], | 
|  | 125 | +   "source": [ | 
|  | 126 | +    "for i in range(len(m_mond)):\n", | 
|  | 127 | +    "    fig, ax =plt.subplots()\n", | 
|  | 128 | +    "    gr.Plot().plot_grc(df_all[i],fit_mond[i],mond,'MOND',galaxy[i],ax)\n", | 
|  | 129 | +    "    gr.Plot().plot_grc(df_all[i],fit_newton[i],newton,'Newton',galaxy[i],ax)\n", | 
|  | 130 | +    "    plt.savefig(galaxy[i]+'_rad_vars_fit.pdf')" | 
|  | 131 | +   ] | 
|  | 132 | +  }, | 
|  | 133 | +  { | 
|  | 134 | +   "cell_type": "code", | 
|  | 135 | +   "execution_count": null, | 
|  | 136 | +   "id": "b4339cbf", | 
|  | 137 | +   "metadata": {}, | 
|  | 138 | +   "outputs": [], | 
|  | 139 | +   "source": [ | 
|  | 140 | +    "fit_mond[1].draw_mnprofile('M0', band=True)" | 
|  | 141 | +   ] | 
|  | 142 | +  }, | 
|  | 143 | +  { | 
|  | 144 | +   "cell_type": "code", | 
|  | 145 | +   "execution_count": null, | 
|  | 146 | +   "id": "d4aa72f7", | 
|  | 147 | +   "metadata": {}, | 
|  | 148 | +   "outputs": [], | 
|  | 149 | +   "source": [ | 
|  | 150 | +    "fit_mond[4].draw_mncontour('rc', 'M0',cl=[.67,.9,.95])" | 
|  | 151 | +   ] | 
|  | 152 | +  } | 
|  | 153 | + ], | 
|  | 154 | + "metadata": { | 
|  | 155 | +  "kernelspec": { | 
|  | 156 | +   "display_name": "Python 3 (ipykernel)", | 
|  | 157 | +   "language": "python", | 
|  | 158 | +   "name": "python3" | 
|  | 159 | +  }, | 
|  | 160 | +  "language_info": { | 
|  | 161 | +   "codemirror_mode": { | 
|  | 162 | +    "name": "ipython", | 
|  | 163 | +    "version": 3 | 
|  | 164 | +   }, | 
|  | 165 | +   "file_extension": ".py", | 
|  | 166 | +   "mimetype": "text/x-python", | 
|  | 167 | +   "name": "python", | 
|  | 168 | +   "nbconvert_exporter": "python", | 
|  | 169 | +   "pygments_lexer": "ipython3", | 
|  | 170 | +   "version": "3.10.6" | 
|  | 171 | +  } | 
|  | 172 | + }, | 
|  | 173 | + "nbformat": 4, | 
|  | 174 | + "nbformat_minor": 5 | 
|  | 175 | +} | 
0 commit comments