Emilys Galactic Gits

Emily's Random Code

Visualization of Astrophysical Data and Fitting Results

SextenFigures.ipynb

This file contains Python code for processing and visualizing astrophysical data, including manipulating photometry data, normalizing densities, creating color-magnitude diagrams, and plotting fitting results to identify best fits among model populations. It includes custom color schemes for visualization, use of matplotlib for plotting, and procedures for density normalization and scatter plot generation. The notebook likely serves as part of an analysis pipeline in astrophysical research, particularly focusing on the study of globular clusters or similar stellar systems.

snippet.python
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import os
import pickle as pkl
from scipy.spatial.distance import cdist
from tqdm import tqdm
import numpy as np
from pathlib import Path
 
from scipy.interpolate import splrep, BSpline
from scipy.signal import savgol_filter
from pysep.atm.utils import load_new_style
from dataclasses import dataclass
from collections import namedtuple
import matplotlib
 
from scipy.signal import savgol_filter
 
from mplEasyAnimate import animation
snippet.python
@dataclass
class colorscheme:
    DartmouthGreen : str = "#00693E"
    ForestGreen : str = "#12312B"
    RichForestGreen : str = "#0D1E1C"
    White : str = "#FFFFFF"
    Black : str = "#000000"
    AutumBrown : str = "#643C20"
    BonfireRed : str = "#9D162E"
    SpringGreen : str = "#C4DD88"
    Violet : str = "#B607DE"
 
    pallet = namedtuple("pallet", "C0 C1")
    primary : namedtuple = pallet(DartmouthGreen, Black)
    seconday : namedtuple = pallet(ForestGreen,AutumBrown)
    tertiery : namedtuple = pallet(RichForestGreen,BonfireRed)
DGcmap =cmap = matplotlib.colors.LinearSegmentedColormap.from_list("", ["#FFFFFF", colorscheme.DartmouthGreen])
snippet.python
PHOTROOT = "/mnt/Astronomy/GraduateSchool/Thesis/GCConsistency/NGC2808/photometry/HUGS/ngc2808/photometry.pkl"
snippet.python
def normalize_density(color, mag, density, n=5000):
    normDensity = np.zeros(shape=color.shape[0])
    for IDx, (c, m, d) in tqdm(enumerate(zip(color, mag, density)), total=len(density)):
        distances = cdist(np.array([[c, m]]), np.array([color, mag]).T)[0]
        closestIDX= np.argpartition(distances, n)
        closestDensity = density[closestIDX[:n]]
        meanNearDensity = np.mean(closestDensity)
        normalizedDensity = d/meanNearDensity
        normDensity[IDx] = normalizedDensity
    return normDensity
def normalize_density_magBin(color, mag, density, binSize=0.1):
    normDensity = np.zeros(shape=color.shape[0])
    for IDx, (c, m, d) in tqdm(enumerate(zip(color, mag, density)), total=len(density)):
        cut = (mag > m-binSize/2) & (mag <= m+binSize/2)
        binDensity = density[cut]
        meanBinDensity = np.mean(binDensity)
        normalizedDensity = d/meanBinDensity
        normDensity[IDx] = normalizedDensity
    return normDensity
snippet.python
densityCache = "/mnt/Astronomy/packages/localTests/fidanka/MC_1_Density.npz"
assert os.path.exists(densityCache), "Density Cache File not Found!"
snippet.python
density = np.load(densityCache)['density']
with open(PHOTROOT, 'rb') as f:
    HUGSPhotometry = pkl.load(f)[1]
color = HUGSPhotometry["F275W"] - HUGSPhotometry["F814W"]
mag = HUGSPhotometry["F814W"]
HUGSPhotometry['density'] = normalize_density_magBin(color, mag, density, binSize=0.3)
100%|████████████████████████████████████████████████████████| 38319/38319 [00:15<00:00, 2418.10it/s]
snippet.python
with plt.style.context(pubStyle):
    fig, axs = plt.subplots(2,3, figsize=(15, 10))
 
    ax1 = axs[0,0]
    ax2 = axs[1,0]
    ax3 = axs[0,1]
    ax4 = axs[1,1]
    ax5 = axs[0,2]
    ax6 = axs[1,2]
 
    f1 = "F275W"
    f2 = "F814W"
    f3 = f2
    color = HUGSPhotometry[f1]-HUGSPhotometry[f2]
    mag = HUGSPhotometry[f3]
    density = HUGSPhotometry["density"]
 
    condF = (color > 1.75) & (color < 5) & (mag < 22) & (mag > 15)
    colorF = color[condF]
    magF = mag[condF]
    densityF = density[condF]
 
    ax1.scatter(colorF, magF, s=1, c=colorscheme.DartmouthGreen)
    ax1.invert_yaxis()
 
    ax2.scatter(colorF, magF, s=1, c=densityF, cmap=DGcmap)
    ax2.invert_yaxis()
 
    condZ = (color > 1.75) & (color < 3.5) & (mag < 20) & (mag > 18.1)
    colorZ = color[condZ]
    magZ = mag[condZ]
    densityZ = density[condZ]
 
    ax3.scatter(colorZ, magZ, s=1, c=colorscheme.DartmouthGreen)
    ax3.set_xticklabels([])
    ax3.invert_yaxis()
 
    # normDensityZ = normalize_density(colorZ, magZ, densityZ)
    ax4.scatter(colorZ, magZ, s=1, c=densityZ, alpha=1, cmap=DGcmap)
    ax4.invert_yaxis()
 
    # IDs = ["A", "B", "C", "D", "E", "F"]
    # Axs = [ax1, ax3, ax5, ax2, ax4, ax6]
    # for ax, l in zip(Axs, IDs):
    #     ax.text(0.9, 0.5, l, transform=ax.transAxes)
 
    plt.subplots_adjust(hspace=0)
 
    ax1.set_ylabel(f"{f3}", fontsize=23)
    ax2.set_ylabel(f"{f3}", fontsize=23)
    # ax3.set_ylabel(f"{f3}", fontsize=23)
    # ax4.set_ylabel(f"{f3}", fontsize=23)
 
    ax2.set_xlabel(f"{f1}-{f2}")
    ax4.set_xlabel(f"{f1}-{f2}")
    ax6.set_xlabel(f"{f1}-{f2}")
 
    condRGB = (color > 3) & (color < 5) & (mag < 17.5) & (mag > 15)
    colorRGB = color[condRGB]
    magRGB = mag[condRGB]
    densityRGB = density[condRGB]
    # normDensity = normalize_density(colorRGB, magRGB, densityRGB, n=10)
    ax5.scatter(colorRGB, magRGB,s=1, c=colorscheme.DartmouthGreen)
    ax5.invert_yaxis()
    ax6.scatter(colorRGB, magRGB,s=1, c=densityRGB, cmap=DGcmap)
    ax6.invert_yaxis()
 
    fig.savefig("Figures/DensityMap.png", dpi=300)

png

snippet.python
with open("ExampleFidOutput.pkl", 'rb') as f:
    fiducial = pkl.load(f)
snippet.python
with plt.style.context(pubStyle):
    fig, ax = plt.subplots(1,1,figsize=(10,7))
    fig.subplots_adjust(hspace=0)
    ax.scatter(color, mag, s=1, c=HUGSPhotometry['density'], cmap=DGcmap)
 
    ax.set_xlim(1.5, 5)
    ax.set_ylim(15.5, 22)
 
    ax.invert_yaxis()
 
    ax.set_xlabel("F275W - F814W", fontsize=23)
    ax.set_ylabel("F814W", fontsize=23)
 
 
    ax.plot(fiducial['Acolor'], fiducial['mag'], linewidth = 5, color=colorscheme.BonfireRed, label="E")
    ax.plot(fiducial['Bcolor'], fiducial['mag'], linewidth = 5, color=colorscheme.Violet, label="A", linestyle='dashed')
 
    ax.legend(fontsize=23, frameon=False)
 
    fig.savefig("Figures/NGC2808Fid.png", dpi=300)
    # fig.savefig("../static/imgs/NGC2808Fid.png", dpi=200)

png

snippet.python
with open("../../../../../../../../Astronomy/GraduateSchool/Thesis/GCConsistency/NGC2808/Analysis/fitting/ReducedResults.denseAlpha.pkl", 'rb') as f:
    fitResults = pkl.load(f)
snippet.python
with plt.style.context(pubStyle):
    for n in tqdm(range(len(fitResults))):
        fig, ax = plt.subplots(1,1,figsize=(10,7))
        ax.plot(fiducial['Acolor'], fiducial['mag'], linewidth = 5, color=colorscheme.BonfireRed, label="E")
        ax.plot(fiducial['Bcolor'], fiducial['mag'], linewidth = 5, color=colorscheme.Violet, label="A", linestyle='dashed')
        ax.set_xlim(1.5, 5)
        ax.set_ylim(15.5, 22)
        # print(fitResults[n][0])
        ax.set_title(f"n: {n}, chi2: {fitResults[n][0]:0.3f}")
        A = fitResults[n][2][0]
        E = fitResults[n][2][1]
 
        Acolor = A['F275W'] - A["F814W"]
        Amag = A["F814W"]
 
        Ecolor = E["F275W"] - E["F814W"]
        Emag = E["F814W"]
        ax.plot(Acolor, Amag)
        ax.plot(Ecolor, Emag)
        ax.invert_yaxis()
 
        plt.savefig(f"Figures/bfr-{n}.png")
        plt.close(fig)
snippet.python
nA = 0 
nE = 4
with plt.style.context(pubStyle):
    A = fitResults[nA][2][0]
    E = fitResults[nE][2][1]
    print(fitResults[nA][0])
    print(fitResults[nE][0])
 
    Acolor = A['WFC3_UVIS_F275W'] - A["WFC3_UVIS_F814W"]
    Amag = A["WFC3_UVIS_F814W"]
 
    Ecolor = E["WFC3_UVIS_F275W"] - E["WFC3_UVIS_F814W"]
    Emag = E["WFC3_UVIS_F814W"]
    fig, ax = plt.subplots(1,1,figsize=(10,7))
 
    ASotedIDX = np.argsort(Amag)
    ax.plot(fiducial['Acolor'], fiducial['mag'], color=colorscheme.BonfireRed, linestyle='dashed', alpha=0.5, label="Fiducial E")
    ax.plot(savgol_filter(Acolor[ASotedIDX], 10, 2), Amag[ASotedIDX], color=colorscheme.BonfireRed, label="Best Fit E")
    ax.plot(fiducial['Bcolor'], fiducial['mag'], color=colorscheme.Violet, linestyle='dashed', alpha=0.5, label="Fiducial A")
    # ax.set_xlim(1.5, 5)
    # ax.set_ylim(18, 21)
    ax.plot(Ecolor, Emag, color=colorscheme.Violet, label="Best Fit A")
    ax.invert_yaxis()
 
    ax.legend(frameon=False)
 
#     print(fitResults[nA][1][0][3].x)
#     print(fitResults[nE][1][1][3].x)
 
#     print(fitResults[nA][1][0][1])
#     print(fitResults[nE][1][1][1])
 
#     print(fitResults[nA][1][0][2])
#     print(fitResults[nE][1][1][2])
 
    ax.set_xlabel("F275W - F814W", fontsize=23)
    ax.set_ylabel("F814W", fontsize=23)
 
    fig.savefig("Figures/BestFit.png", dpi=300)
0.07313541839143663
0.132039504646009

png

snippet.python
A

<div> <style scoped> 

  .dataframe tbody tr th:only-of-type {
      vertical-align: middle;
  }
  .dataframe tbody tr th {
      vertical-align: top;
  }
  .dataframe thead th {
      text-align: right;
  }

</style> <table border=“1” class=“dataframe”

> 

<thead>
  <tr style="text-align: right;">
    <th></th>
    <th>WFC3_UVIS_F200LP</th>
    <th>WFC3_UVIS_F218W</th>
    <th>WFC3_UVIS_F225W</th>
    <th>WFC3_UVIS_F275W</th>
    <th>WFC3_UVIS_F280N</th>
    <th>WFC3_UVIS_F300X</th>
    <th>WFC3_UVIS_F336W</th>
    <th>WFC3_UVIS_F343N</th>
    <th>WFC3_UVIS_F350LP</th>
    <th>WFC3_UVIS_F373N</th>
    <th>...</th>
    <th>WFC3_IR_F127M</th>
    <th>WFC3_IR_F128N</th>
    <th>WFC3_IR_F130N</th>
    <th>WFC3_IR_F132N</th>
    <th>WFC3_IR_F139M</th>
    <th>WFC3_IR_F140W</th>
    <th>WFC3_IR_F153M</th>
    <th>WFC3_IR_F160W</th>
    <th>WFC3_IR_F164N</th>
    <th>WFC3_IR_F167N</th>
  </tr>
</thead>
<tbody>
  <tr>
    <th>0</th>
    <td>25.978353</td>
    <td>33.859351</td>
    <td>33.659541</td>
    <td>32.988767</td>
    <td>30.306583</td>
    <td>30.936987</td>
    <td>29.587342</td>
    <td>29.429675</td>
    <td>25.782097</td>
    <td>29.467213</td>
    <td>...</td>
    <td>23.151844</td>
    <td>23.073539</td>
    <td>23.116001</td>
    <td>23.090822</td>
    <td>22.973418</td>
    <td>22.967202</td>
    <td>22.699599</td>
    <td>22.698120</td>
    <td>22.487449</td>
    <td>22.519156</td>
  </tr>
  <tr>
    <th>1</th>
    <td>25.954552</td>
    <td>33.831236</td>
    <td>33.631871</td>
    <td>32.952047</td>
    <td>30.279008</td>
    <td>30.905088</td>
    <td>29.555683</td>
    <td>29.398252</td>
    <td>25.758344</td>
    <td>29.434568</td>
    <td>...</td>
    <td>23.131067</td>
    <td>23.052768</td>
    <td>23.095286</td>
    <td>23.070104</td>
    <td>22.952688</td>
    <td>22.946525</td>
    <td>22.679091</td>
    <td>22.677667</td>
    <td>22.467302</td>
    <td>22.499040</td>
  </tr>
  <tr>
    <th>2</th>
    <td>25.930738</td>
    <td>33.803047</td>
    <td>33.604014</td>
    <td>32.915213</td>
    <td>30.251354</td>
    <td>30.873142</td>
    <td>29.524005</td>
    <td>29.366815</td>
    <td>25.734579</td>
    <td>29.401911</td>
    <td>...</td>
    <td>23.110263</td>
    <td>23.031970</td>
    <td>23.074543</td>
    <td>23.049356</td>
    <td>22.931943</td>
    <td>22.925829</td>
    <td>22.658564</td>
    <td>22.657194</td>
    <td>22.447123</td>
    <td>22.478889</td>
  </tr>
  <tr>
    <th>3</th>
    <td>25.906920</td>
    <td>33.774848</td>
    <td>33.576146</td>
    <td>32.878354</td>
    <td>30.223691</td>
    <td>30.841181</td>
    <td>29.492311</td>
    <td>29.335363</td>
    <td>25.710810</td>
    <td>29.369246</td>
    <td>...</td>
    <td>23.089460</td>
    <td>23.011173</td>
    <td>23.053802</td>
    <td>23.028610</td>
    <td>22.911199</td>
    <td>22.905135</td>
    <td>22.638039</td>
    <td>22.636723</td>
    <td>22.426947</td>
    <td>22.458742</td>
  </tr>
  <tr>
    <th>4</th>
    <td>25.883110</td>
    <td>33.746654</td>
    <td>33.548278</td>
    <td>32.841488</td>
    <td>30.196032</td>
    <td>30.809220</td>
    <td>29.460617</td>
    <td>29.303914</td>
    <td>25.687049</td>
    <td>29.336588</td>
    <td>...</td>
    <td>23.068670</td>
    <td>22.990389</td>
    <td>23.033073</td>
    <td>23.007876</td>
    <td>22.890467</td>
    <td>22.884452</td>
    <td>22.617526</td>
    <td>22.616263</td>
    <td>22.406783</td>
    <td>22.438606</td>
  </tr>
  <tr>
    <th>...</th>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
    <td>...</td>
  </tr>
  <tr>
    <th>313</th>
    <td>14.949465</td>
    <td>21.859864</td>
    <td>21.407965</td>
    <td>19.522513</td>
    <td>18.613011</td>
    <td>18.281919</td>
    <td>17.055310</td>
    <td>16.928078</td>
    <td>14.766674</td>
    <td>16.485056</td>
    <td>...</td>
    <td>12.671225</td>
    <td>12.594520</td>
    <td>12.638176</td>
    <td>12.607687</td>
    <td>12.484442</td>
    <td>12.482223</td>
    <td>12.211117</td>
    <td>12.225296</td>
    <td>12.052960</td>
    <td>12.080236</td>
  </tr>
  <tr>
    <th>314</th>
    <td>14.916718</td>
    <td>21.845558</td>
    <td>21.411155</td>
    <td>19.530729</td>
    <td>18.592100</td>
    <td>18.274660</td>
    <td>17.044438</td>
    <td>16.916905</td>
    <td>14.733685</td>
    <td>16.470239</td>
    <td>...</td>
    <td>12.631767</td>
    <td>12.554985</td>
    <td>12.598540</td>
    <td>12.567914</td>
    <td>12.444282</td>
    <td>12.441953</td>
    <td>12.169769</td>
    <td>12.183968</td>
    <td>12.010663</td>
    <td>12.037913</td>
  </tr>
  <tr>
    <th>315</th>
    <td>14.884109</td>
    <td>21.831753</td>
    <td>21.415111</td>
    <td>19.539779</td>
    <td>18.571538</td>
    <td>18.267978</td>
    <td>17.034075</td>
    <td>16.906246</td>
    <td>14.700830</td>
    <td>16.455932</td>
    <td>...</td>
    <td>12.592302</td>
    <td>12.515440</td>
    <td>12.558894</td>
    <td>12.528127</td>
    <td>12.404098</td>
    <td>12.401656</td>
    <td>12.128368</td>
    <td>12.142587</td>
    <td>11.968289</td>
    <td>11.995514</td>
  </tr>
  <tr>
    <th>316</th>
    <td>14.851538</td>
    <td>21.818069</td>
    <td>21.419250</td>
    <td>19.549028</td>
    <td>18.551062</td>
    <td>18.261435</td>
    <td>17.023836</td>
    <td>16.895713</td>
    <td>14.668011</td>
    <td>16.441746</td>
    <td>...</td>
    <td>12.552842</td>
    <td>12.475900</td>
    <td>12.519251</td>
    <td>12.488343</td>
    <td>12.363915</td>
    <td>12.361359</td>
    <td>12.086961</td>
    <td>12.101201</td>
    <td>11.925905</td>
    <td>11.953104</td>
  </tr>
  <tr>
    <th>317</th>
    <td>14.818997</td>
    <td>21.804534</td>
    <td>21.423317</td>
    <td>19.559072</td>
    <td>18.530665</td>
    <td>18.255557</td>
    <td>17.014354</td>
    <td>16.885930</td>
    <td>14.635219</td>
    <td>16.428355</td>
    <td>...</td>
    <td>12.513350</td>
    <td>12.436329</td>
    <td>12.479578</td>
    <td>12.448531</td>
    <td>12.323709</td>
    <td>12.321042</td>
    <td>12.045552</td>
    <td>12.059809</td>
    <td>11.883527</td>
    <td>11.910705</td>
  </tr>
</tbody>

</table>

<

p>318 rows × 57 columns</p> </div>

snippet.python
print(type(fitResults))
<class 'list'>
snippet.python
chi2 = list()
Ys = list()
alphas = list()
for pop in fitResults:
    chi2.append(pop[0])
    Ys.append((pop[1][0][1], pop[1][1][1]))
    alphas.append((pop[1][0][2], pop[1][1][2]))
chi2 = np.array(chi2)
Ys = np.array(Ys)
alphas = np.array(alphas)
snippet.python
YsA = Ys[:, 0]
sortedYsA = np.argsort(YsA)
plt.plot(YsA[sortedYsA], chi2[sortedYsA])
plt.xlim(
[<matplotlib.lines.Line2D at 0x7fb4c44179a0>]

png

snippet.python
 
 
 
Ys[:, 0]
array([0.24, 0.3 , 0.24, ..., 0.33, 0.3 , 0.24])
snippet.python
sortedchi2
array([   0,    1,    2, ..., 1079, 1080, 1081])
snippet.python
 
 
 
 
 
    array([0.05120079, 0.06275137, 0.07305412, ..., 0.99586754, 0.99739802,
           0.99955145])