{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "TA3dSyKj1019"
   },
   "source": [
    "# Radio Telescope - Example of Loading Spectrum Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "d0U8GYPu18LH"
   },
   "source": [
    "## Load the Data from a .rad Text File"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "id": "NEmH7StEx4bT"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from scipy import stats\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "SS5gGl5Rx4bU"
   },
   "outputs": [],
   "source": [
    "#%matplotlib notebook"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "1leUFy7oyyz2"
   },
   "source": [
    "Before you run this code, open the file 'Test.rad' using a text editor.  Make sure you are familiar with what is in the file: (more explanation in that file)\n",
    "* rows: different times\n",
    "* columns (after the first few): different frequencies\n",
    "\n",
    "**Task**: Identify which columns you probably will need to load.  For each one, what data format should it be loaded as: string, float, or int?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# (discussion task above)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Review the documentation for np.loadtxt.  Some function keywords that will be useful are dtype, usecols, and comments.  Now, write commands below to load individual columns of the file.  There are partially complete lines to suggest which variables you will need to load."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "id": "tG5EVa4Ax4bV",
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "file=\"Test.rad\"                                              # edit in file name\n",
    "\n",
    "# write these lines yourself (I recommend commenting out all but one and adding one line at a time)\n",
    "time=      # reads time stamp as a string (for now)\n",
    "freq0=        # reads the frequency of the first element in the spectrum\n",
    "freqstep=      # reads the frequency step\n",
    "nspec =           # reads the number of elements in each spectrum\n",
    "\n",
    "# these commands are written for you to set up an argument you can give to usecols\n",
    "endspec=9+nspec[0]  # this is the column position of the last element of each spectrum\n",
    "spec_cols=np.arange(9,endspec,1,dtype=int)            # spectrum is written starting with column 9 and extending for nspec elements\n",
    "\n",
    "# write these lines yourself\n",
    "data=    # data is a two dimensional array - you can feed spec_cols into the usecols argument for loadtxt\n",
    "vlsr=  # VLSR (velocity of local standard of rest) - the column number you need is nspec+1\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "vEDsT0LjKaY2"
   },
   "source": [
    "Freq0 and freqstep should ideally be the same value at every time step - if you change frequency settings, you should save and start a new data file.  Next step is to generate an array containing the list of the frequencies of all the channels (the x-axis for your spectrum plots).\n",
    "\n",
    "You will make plots of your spectrum with frequency on the x-axis and intensity (in \"data\") on the y-axis.  **Task**: Write a command using linspace that defines \"freqval\", the list of frequency values that should go on the x-axis.  freqval should have the same length as spec_cols.  You will need freq0[0], nspec[0], and freqstep[0] to define freqval."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "id": "TRG7_Eawx4bV"
   },
   "outputs": [],
   "source": [
    "# Instructions above - define freqval"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "wiULqJwt1RYH"
   },
   "source": [
    "Add a code cell below.  Check the shape of the 2D data array using data.shape.  Then check the length of the time array and the length of the freqval array.  Out of the rows and columns in data, which direction matches up to frequency, and which to time?  Document your conclusions here."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "WaJKfbsp2B_I"
   },
   "source": [
    "## Plot the Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "kDiyiG9Z_xrR"
   },
   "source": [
    "First plot the data from a single timestamp.  For example, data[7] should contain the spectrum from timestamp #7."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "OFpO8eVGx4bV"
   },
   "outputs": [],
   "source": [
    "# instructions above"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, calculate an average spectrum across all times.  You can do this using np.mean with the axis parameter (look at np.mean's help for more info).  Then, plot the averaged spectrum versus frequency."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "CmFrnWc-x4bW"
   },
   "outputs": [],
   "source": [
    "# instructions above"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "RWVePi8-LIFN"
   },
   "source": [
    "What is the difference between the first and second plot above?  Why are they different? "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "SKxfjjCkLy3b"
   },
   "source": [
    "## Measuring $v_\\textrm{rot}$\n",
    "\n",
    "To measure the galactic rotation speed $v_\\textrm{rot}$, you need to estimate the *maximum* frequency at which there is emission above the baseline level.  Carry out the following steps to get an idea for how you will measure the rotation speed of the galaxy.  You'll eventually collect data (like the above spectrum graph) pointing in many different directions, to make a graph of rotation speed vs. distance from the Galactic Center.\n",
    "\n",
    "1) Estimate this maximum frequency by eye (and start thinking about other more quantitative ways you could measure it later).  What did you get?  This is $f_{obs}$.\n",
    "\n",
    "2) Use $f_{obs}$ to calculate $v_{obs}$.\n",
    "\n",
    "3) Following the explanation in the intro packet, convert this $v_{obs}$ to $v_{los}$.  Since this is the maximum value for $v_{los}$ at which there is emission, we call this the tangential velocity $v_\\textrm{rot}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "1cxC29wtNgOk"
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "colab": {
   "provenance": []
  },
  "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.11.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}