Building ChRMDir excel file from csv.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "10920f38",
   "metadata": {},
   "source": [
    "# Building  Pmag DIR excel file from Paldir's *.csv files"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4f4da139",
   "metadata": {},
   "source": [
    "### Import pmag modules and function definitions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "7b6d38a7",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pmagpy.pmag as pmag\n",
    "import pmagpy.pmagplotlib as pmagplotlib\n",
    "import pmagpy.ipmag as ipmag\n",
    "import matplotlib.pyplot as plt # our plotting buddy\n",
    "from pmagpy import convert_2_magic as convert\n",
    "import numpy as np # the fabulous NumPy package\n",
    "import pandas as pd # and of course Pandas\n",
    "has_basemap, Basemap = pmag.import_basemap()\n",
    "has_cartopy, Cartopy = pmag.import_cartopy()\n",
    "from IPython.display import Image\n",
    "%matplotlib inline\n",
    "\n",
    "#define a xlookup function\n",
    "def xlookup(lookup_value, lookup_array, return_array, if_not_found:str = ''):\n",
    "    match_value = return_array.loc[lookup_array == lookup_value]\n",
    "    if match_value.empty:\n",
    "        return f'\"{lookup_value}\" not found!' if if_not_found == '' else if_not_found\n",
    "\n",
    "    else:\n",
    "        return match_value.tolist()[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e53eae1a",
   "metadata": {},
   "source": [
    "#### Set geographic latitude and longitude of the site (in degrees)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "244809d4",
   "metadata": {},
   "outputs": [],
   "source": [
    "lat = 37.72 # ENTER site latitude\n",
    "lon = -2.47 # ENTER site longitude"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fc7a50d5",
   "metadata": {},
   "source": [
    "### Selecting input files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "f1feac09",
   "metadata": {},
   "outputs": [],
   "source": [
    "path = ''          #Give a path to folder from working directory\n",
    "site_name = \"\"       #Give a site name\n",
    "\n",
    "file1 = path+\"\"     #Give csv file names (directions picked with 3gpaldir program). EX: BH1.CSV\n",
    "file2 = path+\"\"\n",
    "thickness_filename = '.xlsx' # ex: 'borreda_level.xlsx'\n",
    "bedding_filename = '' # ex: borreda_bedding.xlsx\n",
    "\n",
    "thickness_file = path+thickness_filename #give stratigraphic position in meters\n",
    "bedding_file = path+bedding_filename #give bedding orientation (dip direction and dip)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2c7dbe9a",
   "metadata": {},
   "source": [
    "#### Build dataframe from csv files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "7cfdc056",
   "metadata": {},
   "outputs": [],
   "source": [
    "header= [\"sample\",\"level (m)\",\"dec\",\"inc\",\"Int(E-6A/m)\",\"error\",\"Q\",\"Temp\",\"comments\"]\n",
    "df_file1 = pd.read_csv(file1, names = header)\n",
    "\n",
    "if file2:\n",
    "    df_file2 = pd.read_csv(file2, names = header)\n",
    "    df_vectors = df_file1.append(df_file2)\n",
    "else:\n",
    "    df_vectors = df_file1\n",
    "    \n",
    "df_vectors"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "39c86add",
   "metadata": {},
   "source": [
    "#### Read sample level (stratigraphic thickness)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "b201926d",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Set the number of characters to read from sample name (default n = 5): 4\n"
     ]
    }
   ],
   "source": [
    "sample_characters = str(input('Set the number of characters to read from sample name (default n = 5): '))\n",
    "if not sample_characters:\n",
    "    sample_characters = 5\n",
    "else:\n",
    "    sample_characters = int(sample_characters)\n",
    "    \n",
    "#adjust site-names to selected sample_characters number \n",
    "df_vectors['newsamplename'] = df_vectors['sample'].str[:sample_characters] #taking the first 4 characters from left side.\n",
    "\n",
    "#reading data from excel with column labels: 'sample' and 'level'\n",
    "#path+\"borreda_level.xlsx\"\n",
    "if thickness_filename:\n",
    "    xls = pd.ExcelFile(thickness_file)\n",
    "    samplelevel = xls.parse(site_name).round(1)\n",
    "    cols = samplelevel.columns.tolist()\n",
    "    samplelevel = samplelevel.rename(columns = {cols[0]:'sample', cols[1]:'level (m)'})\n",
    "    # take the first 4 characters of sample name\n",
    "    samplelevel['newsamplename'] = samplelevel['sample'].str[:sample_characters] #taking the first 4 characters from left side.\n",
    "    #samplelevel['sample'] =samplelevel['newsamplename']\n",
    "    #samplelevel = samplelevel.drop(columns=['newsamplename'])\n",
    "    #samplelevel.dropna(subset=['sample'], inplace=True)\n",
    "    #samplelevel = samplelevel.reset_index(drop = True)  \n",
    "    # ADD stratigraphic level to df_vectors dataframe\n",
    "    df_vectors['level (m)'] = df_vectors['newsamplename'].apply(xlookup, args = (samplelevel['newsamplename'], samplelevel['level (m)']))\n",
    "\n",
    "df_vectors"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "20f70a05",
   "metadata": {},
   "source": [
    "#### Read bedding orientation file (uncomment ''' ''' if any)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "5295bbd6",
   "metadata": {},
   "outputs": [],
   "source": [
    "'''\n",
    "#reading data from excel spreadsheets and importing into a single dataframe\n",
    "if bedding_filename:\n",
    "    xls = pd.ExcelFile(bedding_file)\n",
    "    dipdir_dip = xls.parse(0).round(1)\n",
    "\n",
    "    cols = dipdir_dip.columns.tolist()\n",
    "    dipdir_dip = dipdir_dip.rename(columns = {cols[0]:'sample', cols[1]:'So_dipdir', cols[2]:'So_dip'})\n",
    "    dipdir_dip['newsamplename'] = dipdir_dip['sample'].str[:sample_characters] #taking the first 4 characters from left side.\n",
    "    # ADD So_dipdir to df_vectors dataframe\n",
    "    df_vectors['So_dipdir'] = df_vectors['newsamplename'].apply(xlookup, args = (dipdir_dip['newsamplename'], dipdir_dip['So_dipdir']))\n",
    "    # ADD So_dip to df_vectors dataframe\n",
    "    df_vectors['So_dip'] = df_vectors['newsamplename'].apply(xlookup, args = (dipdir_dip['newsamplename'], dipdir_dip['So_dip']))\n",
    "\n",
    "    dipdir_dip.head()\n",
    "'''"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "da47abf1",
   "metadata": {},
   "outputs": [],
   "source": [
    "df_vectors = df_vectors.drop(columns=['newsamplename'])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ccb6830f",
   "metadata": {},
   "source": [
    "#### Format dataframe df_vectors (only if vectors are in geographic and stratigraphic coordinates)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "42e2aabb",
   "metadata": {},
   "outputs": [],
   "source": [
    "# SKIP THIS BLOCK IF ONLY ONE COORDINATE SYSTEM\n",
    "\n",
    "df_vectors['Dg'] = ''\n",
    "df_vectors['Ig'] = ''\n",
    "df_vectors['Ds'] = ''\n",
    "df_vectors['Is'] = ''\n",
    "comment = df_vectors['comments']\n",
    "intensity = df_vectors['Int(E-6A/m)']\n",
    "quality = df_vectors['Q']\n",
    "sample = df_vectors['sample']\n",
    "\n",
    "#loop to fill (Dg, Ig) and (Ds, Is) columns\n",
    "\n",
    "for value in range(len(df_vectors)):\n",
    "    if 'notc' in str(comment[value]):\n",
    "        df_vectors.loc[value,'Dg'] = df_vectors['dec'][value]\n",
    "        df_vectors.loc[value,'Ig'] = df_vectors['inc'][value]\n",
    "        for value2 in range(len(df_vectors)):\n",
    "            if ' tc' in str(comment[value2]):\n",
    "                if sample[value] == sample[value2] and quality[value] == quality[value2]:\n",
    "                    df_vectors.loc[value,'Ds'] = df_vectors['dec'][value2]\n",
    "                    df_vectors.loc[value,'Is'] = df_vectors['inc'][value2]\n",
    "    elif 'tc' in str(comment[value]):\n",
    "        df_vectors.loc[value,'Ds'] = df_vectors['dec'][value]\n",
    "        df_vectors.loc[value,'Is'] = df_vectors['inc'][value]\n",
    "\n",
    "df_vectors['Dg'].replace('', np.nan, inplace=True)      \n",
    "df_vectors['Ig'].replace('', np.nan, inplace=True)\n",
    "df_vectors.dropna(subset=['sample'], inplace=True)\n",
    "df_vectors = df_vectors.reset_index(drop = True)\n",
    "df_vectors = df_vectors.replace(r'^\\s*$', np.nan, regex=True)   # replace empty space with NaN\n",
    "df_vectors = df_vectors.dropna(subset=['Dg'])                   # delete rows with NaN in column Dg\n",
    "df_vectors =df_vectors.dropna(axis='columns', how= 'all')       # delete empty columns\n",
    "df_vectors[\"comments\"] = df_vectors[\"comments\"].str.replace(\"\\snotc\", \"\")  # delete string \"notc\" from columns comments\n",
    "\n",
    "#re-arrange columns of dataframe\n",
    "cols = df_vectors.columns.tolist()\n",
    "cols = cols[0:2] + cols[9:] + cols[4:9]\n",
    "df_vectors = df_vectors[cols]                  #reordering columns\n",
    "df_vectors = df_vectors.round(2)\n",
    "df_vectors = df_vectors.set_index('sample')        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e012c7df",
   "metadata": {},
   "outputs": [],
   "source": [
    "'''\n",
    "# EXECUTE THIS LINES IF THE BLOCK ABOVE WAS NOT EXECUTED\n",
    "df_vectors = df_vectors.drop(columns=['newsamplename'])  #drop unnecessary columns\n",
    "df_vectors = df_vectors.rename(columns = {'dec': 'Ds', 'inc': 'Is'})\n",
    "df_vectors\n",
    "'''"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b325c948",
   "metadata": {},
   "source": [
    "### VGP calculation from Ds Is and site coordinates"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "ea081d23",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Building np.array\n",
    "chrm_di =np.array([df_vectors[\"Ds\"], df_vectors[\"Is\"]]).transpose()\n",
    "Lats = np.full((len(chrm_di)),lat)\n",
    "Lons = np.full((len(chrm_di)),lon)\n",
    "a95=np.zeros(len(chrm_di))\n",
    "DecIncLatLons = np.column_stack((chrm_di,a95,Lats,Lons))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "5f5f946e",
   "metadata": {},
   "outputs": [],
   "source": [
    "#help(pmag.dia_vgp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "236cb1b3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>level (m)</th>\n",
       "      <th>Dg</th>\n",
       "      <th>Ig</th>\n",
       "      <th>Ds</th>\n",
       "      <th>Is</th>\n",
       "      <th>Int(E-6A/m)</th>\n",
       "      <th>error</th>\n",
       "      <th>Q</th>\n",
       "      <th>Temp</th>\n",
       "      <th>comments</th>\n",
       "      <th>VGP lat</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>sample</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>BH02A</th>\n",
       "      <td>2.2</td>\n",
       "      <td>192.13</td>\n",
       "      <td>-11.77</td>\n",
       "      <td>184.94</td>\n",
       "      <td>-23.89</td>\n",
       "      <td>244.33</td>\n",
       "      <td>6.74</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-240</td>\n",
       "      <td>-64.4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH03-2</th>\n",
       "      <td>4.1</td>\n",
       "      <td>10.89</td>\n",
       "      <td>32.83</td>\n",
       "      <td>352.65</td>\n",
       "      <td>42.35</td>\n",
       "      <td>56.35</td>\n",
       "      <td>2.87</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-240</td>\n",
       "      <td>75.4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH04-2A</th>\n",
       "      <td>5.1</td>\n",
       "      <td>23.64</td>\n",
       "      <td>52.49</td>\n",
       "      <td>345.92</td>\n",
       "      <td>63.70</td>\n",
       "      <td>88.37</td>\n",
       "      <td>2.19</td>\n",
       "      <td>2</td>\n",
       "      <td>5</td>\n",
       "      <td>240-280</td>\n",
       "      <td>77.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH06A A</th>\n",
       "      <td>8.8</td>\n",
       "      <td>216.93</td>\n",
       "      <td>-56.61</td>\n",
       "      <td>172.23</td>\n",
       "      <td>-72.11</td>\n",
       "      <td>731.09</td>\n",
       "      <td>8.67</td>\n",
       "      <td>1</td>\n",
       "      <td>3</td>\n",
       "      <td>240-410</td>\n",
       "      <td>-69.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH05 2A</th>\n",
       "      <td>6.0</td>\n",
       "      <td>13.05</td>\n",
       "      <td>41.12</td>\n",
       "      <td>348.50</td>\n",
       "      <td>50.33</td>\n",
       "      <td>249.97</td>\n",
       "      <td>2.13</td>\n",
       "      <td>2</td>\n",
       "      <td>4</td>\n",
       "      <td>200-280</td>\n",
       "      <td>78.4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH07-1</th>\n",
       "      <td>11.4</td>\n",
       "      <td>156.97</td>\n",
       "      <td>-32.81</td>\n",
       "      <td>142.21</td>\n",
       "      <td>-28.18</td>\n",
       "      <td>510.05</td>\n",
       "      <td>6.62</td>\n",
       "      <td>2</td>\n",
       "      <td>4</td>\n",
       "      <td>320-380</td>\n",
       "      <td>-49.7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH08-2</th>\n",
       "      <td>12.7</td>\n",
       "      <td>261.04</td>\n",
       "      <td>-16.42</td>\n",
       "      <td>264.05</td>\n",
       "      <td>-40.85</td>\n",
       "      <td>443.41</td>\n",
       "      <td>5.21</td>\n",
       "      <td>2</td>\n",
       "      <td>4</td>\n",
       "      <td>240-380</td>\n",
       "      <td>-18.5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH09 2</th>\n",
       "      <td>14.8</td>\n",
       "      <td>351.47</td>\n",
       "      <td>52.93</td>\n",
       "      <td>318.85</td>\n",
       "      <td>50.70</td>\n",
       "      <td>4267.00</td>\n",
       "      <td>2.99</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>320-560</td>\n",
       "      <td>55.8</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH10A</th>\n",
       "      <td>16.3</td>\n",
       "      <td>238.29</td>\n",
       "      <td>-56.85</td>\n",
       "      <td>209.83</td>\n",
       "      <td>-80.10</td>\n",
       "      <td>984.95</td>\n",
       "      <td>6.67</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>200-350</td>\n",
       "      <td>-53.5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH11-1A</th>\n",
       "      <td>17.8</td>\n",
       "      <td>181.50</td>\n",
       "      <td>-32.54</td>\n",
       "      <td>163.85</td>\n",
       "      <td>-38.18</td>\n",
       "      <td>112.98</td>\n",
       "      <td>4.76</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-240</td>\n",
       "      <td>-68.6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH11-2</th>\n",
       "      <td>17.8</td>\n",
       "      <td>269.69</td>\n",
       "      <td>-24.66</td>\n",
       "      <td>277.07</td>\n",
       "      <td>-47.71</td>\n",
       "      <td>121.75</td>\n",
       "      <td>5.01</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-240</td>\n",
       "      <td>-12.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH13A</th>\n",
       "      <td>22.3</td>\n",
       "      <td>183.71</td>\n",
       "      <td>-1.70</td>\n",
       "      <td>181.02</td>\n",
       "      <td>-11.35</td>\n",
       "      <td>40.33</td>\n",
       "      <td>8.70</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>240-280</td>\n",
       "      <td>-58.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH15</th>\n",
       "      <td>25.6</td>\n",
       "      <td>172.09</td>\n",
       "      <td>-49.53</td>\n",
       "      <td>160.16</td>\n",
       "      <td>-50.61</td>\n",
       "      <td>45.66</td>\n",
       "      <td>5.89</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-280</td>\n",
       "      <td>-72.5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH16A</th>\n",
       "      <td>28.6</td>\n",
       "      <td>217.05</td>\n",
       "      <td>-17.97</td>\n",
       "      <td>214.77</td>\n",
       "      <td>-26.26</td>\n",
       "      <td>252.48</td>\n",
       "      <td>2.73</td>\n",
       "      <td>2</td>\n",
       "      <td>4</td>\n",
       "      <td>200-280</td>\n",
       "      <td>-51.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH17</th>\n",
       "      <td>29.5</td>\n",
       "      <td>191.74</td>\n",
       "      <td>-27.49</td>\n",
       "      <td>186.70</td>\n",
       "      <td>-32.38</td>\n",
       "      <td>245.79</td>\n",
       "      <td>5.57</td>\n",
       "      <td>2</td>\n",
       "      <td>5</td>\n",
       "      <td>200-280</td>\n",
       "      <td>-69.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH18</th>\n",
       "      <td>30.4</td>\n",
       "      <td>159.76</td>\n",
       "      <td>-50.24</td>\n",
       "      <td>147.99</td>\n",
       "      <td>-49.17</td>\n",
       "      <td>115.79</td>\n",
       "      <td>3.71</td>\n",
       "      <td>2</td>\n",
       "      <td>4</td>\n",
       "      <td>200-280</td>\n",
       "      <td>-62.5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH19</th>\n",
       "      <td>35.8</td>\n",
       "      <td>168.56</td>\n",
       "      <td>-58.19</td>\n",
       "      <td>152.32</td>\n",
       "      <td>-58.26</td>\n",
       "      <td>152.62</td>\n",
       "      <td>5.58</td>\n",
       "      <td>2</td>\n",
       "      <td>5</td>\n",
       "      <td>200-280</td>\n",
       "      <td>-68.3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH20B</th>\n",
       "      <td>37.4</td>\n",
       "      <td>161.12</td>\n",
       "      <td>-41.11</td>\n",
       "      <td>152.46</td>\n",
       "      <td>-40.55</td>\n",
       "      <td>1265.21</td>\n",
       "      <td>3.30</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>240-380</td>\n",
       "      <td>-62.3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH21A</th>\n",
       "      <td>38.6</td>\n",
       "      <td>164.37</td>\n",
       "      <td>13.68</td>\n",
       "      <td>166.71</td>\n",
       "      <td>12.72</td>\n",
       "      <td>152.66</td>\n",
       "      <td>2.10</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>240-280</td>\n",
       "      <td>-44.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH22A</th>\n",
       "      <td>40.2</td>\n",
       "      <td>6.32</td>\n",
       "      <td>18.96</td>\n",
       "      <td>2.81</td>\n",
       "      <td>23.12</td>\n",
       "      <td>31.24</td>\n",
       "      <td>11.87</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-240</td>\n",
       "      <td>64.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH23</th>\n",
       "      <td>46.2</td>\n",
       "      <td>353.87</td>\n",
       "      <td>-5.30</td>\n",
       "      <td>354.56</td>\n",
       "      <td>-2.83</td>\n",
       "      <td>67.51</td>\n",
       "      <td>14.93</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-280</td>\n",
       "      <td>50.5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH24B</th>\n",
       "      <td>47.7</td>\n",
       "      <td>341.22</td>\n",
       "      <td>53.69</td>\n",
       "      <td>327.84</td>\n",
       "      <td>52.73</td>\n",
       "      <td>67.27</td>\n",
       "      <td>7.27</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>240-280</td>\n",
       "      <td>63.6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH25B</th>\n",
       "      <td>48.7</td>\n",
       "      <td>190.27</td>\n",
       "      <td>69.26</td>\n",
       "      <td>207.87</td>\n",
       "      <td>62.87</td>\n",
       "      <td>38.45</td>\n",
       "      <td>13.20</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>240-280</td>\n",
       "      <td>-4.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH27B</th>\n",
       "      <td>50.6</td>\n",
       "      <td>313.42</td>\n",
       "      <td>57.08</td>\n",
       "      <td>301.60</td>\n",
       "      <td>51.67</td>\n",
       "      <td>77.10</td>\n",
       "      <td>6.75</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-280</td>\n",
       "      <td>42.6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH28A</th>\n",
       "      <td>51.4</td>\n",
       "      <td>350.95</td>\n",
       "      <td>63.03</td>\n",
       "      <td>341.05</td>\n",
       "      <td>63.55</td>\n",
       "      <td>31.88</td>\n",
       "      <td>4.02</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>240-280</td>\n",
       "      <td>74.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH29A</th>\n",
       "      <td>52.4</td>\n",
       "      <td>353.40</td>\n",
       "      <td>63.97</td>\n",
       "      <td>343.08</td>\n",
       "      <td>64.69</td>\n",
       "      <td>103.34</td>\n",
       "      <td>4.49</td>\n",
       "      <td>2</td>\n",
       "      <td>5</td>\n",
       "      <td>200-280</td>\n",
       "      <td>74.7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH31A</th>\n",
       "      <td>53.8</td>\n",
       "      <td>32.29</td>\n",
       "      <td>21.73</td>\n",
       "      <td>30.93</td>\n",
       "      <td>25.65</td>\n",
       "      <td>46.90</td>\n",
       "      <td>9.38</td>\n",
       "      <td>2</td>\n",
       "      <td>5</td>\n",
       "      <td>240-280</td>\n",
       "      <td>53.4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH32A</th>\n",
       "      <td>54.7</td>\n",
       "      <td>61.72</td>\n",
       "      <td>49.05</td>\n",
       "      <td>60.76</td>\n",
       "      <td>53.99</td>\n",
       "      <td>36.87</td>\n",
       "      <td>2.61</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>240-280</td>\n",
       "      <td>41.7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH38A</th>\n",
       "      <td>58.2</td>\n",
       "      <td>354.68</td>\n",
       "      <td>13.00</td>\n",
       "      <td>353.51</td>\n",
       "      <td>14.21</td>\n",
       "      <td>48.26</td>\n",
       "      <td>2.12</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-240</td>\n",
       "      <td>58.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH43</th>\n",
       "      <td>62.8</td>\n",
       "      <td>327.88</td>\n",
       "      <td>74.17</td>\n",
       "      <td>312.15</td>\n",
       "      <td>72.45</td>\n",
       "      <td>229.87</td>\n",
       "      <td>1.95</td>\n",
       "      <td>3</td>\n",
       "      <td>5</td>\n",
       "      <td>200-240</td>\n",
       "      <td>53.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>BH44</th>\n",
       "      <td>73.0</td>\n",
       "      <td>359.33</td>\n",
       "      <td>43.03</td>\n",
       "      <td>354.75</td>\n",
       "      <td>44.50</td>\n",
       "      <td>1169.18</td>\n",
       "      <td>10.92</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>240-530</td>\n",
       "      <td>77.6</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "            level (m)      Dg     Ig      Ds     Is  Int(E-6A/m)  error  Q  \\\n",
       "sample                                                                       \n",
       "BH02A             2.2  192.13 -11.77  184.94 -23.89       244.33   6.74  3   \n",
       "BH03-2            4.1   10.89  32.83  352.65  42.35        56.35   2.87  3   \n",
       "BH04-2A           5.1   23.64  52.49  345.92  63.70        88.37   2.19  2   \n",
       "BH06A A           8.8  216.93 -56.61  172.23 -72.11       731.09   8.67  1   \n",
       "BH05 2A           6.0   13.05  41.12  348.50  50.33       249.97   2.13  2   \n",
       "BH07-1           11.4  156.97 -32.81  142.21 -28.18       510.05   6.62  2   \n",
       "BH08-2           12.7  261.04 -16.42  264.05 -40.85       443.41   5.21  2   \n",
       "BH09 2           14.8  351.47  52.93  318.85  50.70      4267.00   2.99  1   \n",
       "BH10A            16.3  238.29 -56.85  209.83 -80.10       984.95   6.67  1   \n",
       "BH11-1A          17.8  181.50 -32.54  163.85 -38.18       112.98   4.76  3   \n",
       "BH11-2           17.8  269.69 -24.66  277.07 -47.71       121.75   5.01  3   \n",
       "BH13A            22.3  183.71  -1.70  181.02 -11.35        40.33   8.70  3   \n",
       "BH15             25.6  172.09 -49.53  160.16 -50.61        45.66   5.89  3   \n",
       "BH16A            28.6  217.05 -17.97  214.77 -26.26       252.48   2.73  2   \n",
       "BH17             29.5  191.74 -27.49  186.70 -32.38       245.79   5.57  2   \n",
       "BH18             30.4  159.76 -50.24  147.99 -49.17       115.79   3.71  2   \n",
       "BH19             35.8  168.56 -58.19  152.32 -58.26       152.62   5.58  2   \n",
       "BH20B            37.4  161.12 -41.11  152.46 -40.55      1265.21   3.30  1   \n",
       "BH21A            38.6  164.37  13.68  166.71  12.72       152.66   2.10  3   \n",
       "BH22A            40.2    6.32  18.96    2.81  23.12        31.24  11.87  3   \n",
       "BH23             46.2  353.87  -5.30  354.56  -2.83        67.51  14.93  3   \n",
       "BH24B            47.7  341.22  53.69  327.84  52.73        67.27   7.27  3   \n",
       "BH25B            48.7  190.27  69.26  207.87  62.87        38.45  13.20  3   \n",
       "BH27B            50.6  313.42  57.08  301.60  51.67        77.10   6.75  3   \n",
       "BH28A            51.4  350.95  63.03  341.05  63.55        31.88   4.02  3   \n",
       "BH29A            52.4  353.40  63.97  343.08  64.69       103.34   4.49  2   \n",
       "BH31A            53.8   32.29  21.73   30.93  25.65        46.90   9.38  2   \n",
       "BH32A            54.7   61.72  49.05   60.76  53.99        36.87   2.61  3   \n",
       "BH38A            58.2  354.68  13.00  353.51  14.21        48.26   2.12  3   \n",
       "BH43             62.8  327.88  74.17  312.15  72.45       229.87   1.95  3   \n",
       "BH44             73.0  359.33  43.03  354.75  44.50      1169.18  10.92  1   \n",
       "\n",
       "            Temp    comments  VGP lat  \n",
       "sample                                 \n",
       "BH02A          5     200-240    -64.4  \n",
       "BH03-2         5     200-240     75.4  \n",
       "BH04-2A        5     240-280     77.0  \n",
       "BH06A A        3     240-410    -69.9  \n",
       "BH05 2A        4     200-280     78.4  \n",
       "BH07-1         4     320-380    -49.7  \n",
       "BH08-2         4     240-380    -18.5  \n",
       "BH09 2         2     320-560     55.8  \n",
       "BH10A          2     200-350    -53.5  \n",
       "BH11-1A        5     200-240    -68.6  \n",
       "BH11-2         5     200-240    -12.1  \n",
       "BH13A          5     240-280    -58.0  \n",
       "BH15           5     200-280    -72.5  \n",
       "BH16A          4     200-280    -51.0  \n",
       "BH17           5     200-280    -69.0  \n",
       "BH18           4     200-280    -62.5  \n",
       "BH19           5     200-280    -68.3  \n",
       "BH20B          2     240-380    -62.3  \n",
       "BH21A          5     240-280    -44.1  \n",
       "BH22A          5     200-240     64.2  \n",
       "BH23           5     200-280     50.5  \n",
       "BH24B          5     240-280     63.6  \n",
       "BH25B          5     240-280     -4.2  \n",
       "BH27B          5     200-280     42.6  \n",
       "BH28A          5     240-280     74.0  \n",
       "BH29A          5     200-280     74.7  \n",
       "BH31A          5     240-280     53.4  \n",
       "BH32A          5     240-280     41.7  \n",
       "BH38A          5     200-240     58.9  \n",
       "BH43           5     200-240     53.2  \n",
       "BH44           2     240-530     77.6  "
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Run dia_vgp pmag function to calculate VGP latitude\n",
    "vgps=np.array(pmag.dia_vgp(DecIncLatLons)) # get a tuple with lat,lon,dp,dm, convert to array\n",
    "# print(vgps[1:2].transpose()) #  print out the vgps\n",
    "df_vectors['VGP lat'] = vgps[1:2].transpose().round(1)\n",
    "df_vectors"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f4d48727",
   "metadata": {},
   "source": [
    "### Export  DIR file to excel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "3de28441",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "export Huescar_Dir.xlsx (overwrites)? (y/n):y\n",
      "file saved to  Huescar _Dir.xlsx\n"
     ]
    }
   ],
   "source": [
    "outputfile = path+site_name+\"_Dir.xlsx\"\n",
    "export_file = ''\n",
    "while export_file != \"y\" and export_file != \"n\":\n",
    "    export_file = str(input(\"export \"+outputfile+\" (overwrites)? (y/n):\"))\n",
    "    if export_file == 'y':\n",
    "        df_vectors.to_excel(outputfile)\n",
    "        print(\"file saved to \", path+site_name,\"_Dir.xlsx\")\n",
    "    if export_file == \"n\":\n",
    "        df_vectors.to_excel(path+'output_Dir.xlsx')\n",
    "        print(\"file saved to \", path+\"output_Dir.xlsx\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9f4c6c14",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}