# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 3
#     language: python
#     name: python3
# ---

import sys
print(sys.version)

"""
Created on January 30 2019
@author: <NAME>
@contact: <EMAIL>
@web: www.ncaplar.com
"""

# + language="javascript"
# try {
#   require(['base/js/utils'], function (utils) {
#     utils.load_extension('code_prettify/code_prettify');
#     utils.load_extension('collapsible_headings/main'); 
#     utils.load_extension('codefolding/edit'); 
#     utils.load_extension('codefolding/main'); 
#     utils.load_extension('execute_time/ExecuteTime');   
#     utils.load_extension('toc2/main'); 
#   });
# }
# catch (err) {
#   console.log('toc2 load error:', err);
# }
# -

# make notebook nice and wide to fill the entire screen
from IPython.core.display import display, HTML
display(HTML("<style>.container { width:100% !important; }</style>"))

# +
import numpy as np
np.set_printoptions(suppress=True)
np.seterr(divide='ignore', invalid='ignore')

# astropy
from astropy.io import *
from astropy.io import fits

import pickle

#matplotlib
import matplotlib
import matplotlib.pyplot as plt
from matplotlib.colors import LogNorm
matplotlib.rcParams.update({'font.size': 18})
# %config InlineBackend.rc = {}
# %matplotlib inline
# -

# # Original detector map and observed data

# +
detectorMapOriginal=fits.open('/Users/nevencaplar/Documents/PFS/DetectorMap/detectorMap-sim-1-r.fits')
detectorMap5833=fits.open('/Users/nevencaplar/Documents/PFS/DetectorMap/pfsDetectorMap-005833-r1.fits')


hdu = detectorMapOriginal["FIBERID"]
fiberIds = hdu.data
fiberIds = fiberIds.astype(np.int32) 
fiberIds.shape

list_of_illuminated_slit_fibers_in_January_24_data=[32,63,111,192,223,255,289,308,342,401,418,464,518,525,587,620]

# +
list_of_illuminated_science_fibers_in_January_24_data=[]
for i in list_of_illuminated_slit_fibers_in_January_24_data:
    list_of_illuminated_science_fibers_in_January_24_data.append(np.where(fiberIds==i)[0][0])
    
list_of_illuminated_science_fibers_in_January_24_data

# +
with open('/Users/nevencaplar/Documents/PFS/Data_Nov_14/Dataframes/finalHgAr_expanded.pkl', 'rb') as f:
    finalHgAr=pickle.load(f)
    
# I know this was the brightest spot in the old data
finalHgAr[finalHgAr["old_index"]==25]
# -

detectorMap5833[3].data[599][1]

plt.figure(figsize=(8,8))
plt.scatter(detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list[0][:,1], detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list[0][:,2], s=15, c='red', cmap='seismic', alpha=1,edgecolor='black')
plt.ylim(0,4176)
plt.xlim(0,4096)
plt.xlabel('x position on chip')
plt.ylabel('y position on chip')
plt.title('Original DetectorMap' )

detectorMapOriginal[3].data

detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_x=[]
detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_y=[]
for j in list_of_illuminated_science_fibers_in_January_24_data:
    print(j)
    detectorMapOriginal_single_fiber_x=detectorMapOriginal[3].data[detectorMapOriginal[3].data['index']==j]    
    detectorMapOriginal_single_fiber_y=detectorMapOriginal[4].data[detectorMapOriginal[4].data['index']==j]
    
    detectorMapOriginal_single_fiber_as_list_x=[]
    detectorMapOriginal_single_fiber_as_list_y=[]
    # because I do not know how to work with numpy.record
    for i in range(len(detectorMapOriginal_single_fiber_y)):
        detectorMapOriginal_single_fiber_as_list_x.append(list(detectorMapOriginal_single_fiber_x[i]))
        detectorMapOriginal_single_fiber_as_list_y.append(list(detectorMapOriginal_single_fiber_y[i]))
        
    detectorMapOriginal_single_fiber_as_array_x=np.array(detectorMapOriginal_single_fiber_as_list_x)
    detectorMapOriginal_single_fiber_as_array_y=np.array(detectorMapOriginal_single_fiber_as_list_y)
    detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_x.append(detectorMapOriginal_single_fiber_as_array_x)
    detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_y.append(detectorMapOriginal_single_fiber_as_array_y)


# +
def find_nearest(array, value):
    array = np.asarray(array)
    idx = (np.abs(array - value)).argmin()
    return array[idx]

x_y_values_for_912=[]
for i in range(16):
    value_nearest=find_nearest(detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_y[i][:,2],912.58)
    y_value_for_912_single_fiber=detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_y[i][detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_y[i][:,2]==value_nearest]
    
    value_nearest=find_nearest(detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_x[i][:,1],y_value_for_912_single_fiber[0][1])
    x_value_for_912_single_fiber=detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_x[i][detectorMapOriginal_for_illuminated_fiber_Jan_24_2019_as_list_x[i][:,1]==value_nearest]
    
    
    x_y_values_for_912.append([x_value_for_912_single_fiber,y_value_for_912_single_fiber])
# -

x_y_values_for_912

# +
x_y_values_for_912_concise=[]
for i in range(16):
    x_y_values_for_912_concise.append([x_y_values_for_912[i][0][0][2],x_y_values_for_912[i][1][0][1]])
    
x_y_values_for_912_concise=np.array(x_y_values_for_912_concise)


list_of_x_coordinates=[272,486,679,726,1329,1606,1712,2086,2296,2412,2618,2998,3192,3595,4132,4325]
list_of_y_coordinates=[3414,3408,3391,3391,3384,3381,3378,3378,3378,3378,3380,3387,3387,3407,3433,3420]

plt.figure(figsize=(8,8))
plt.scatter(x_y_values_for_912_concise[:,0], x_y_values_for_912_concise[:,1], s=15, c='red', cmap='seismic', alpha=1,edgecolor='black')
plt.scatter(list_of_x_coordinates, list_of_y_coordinates, s=15, c='blue', cmap='seismic', alpha=1,edgecolor='black')
plt.ylim(0,4176)
plt.xlim(0,4406)
plt.xlabel('x position on chip')
plt.ylabel('y position on chip')
plt.title('Original DetectorMap' )
# -

plt.figure(figsize=(12,8))
plt.plot([3415,3405,3391,3391,3384,3383,3379,3378,3378,3378,3380,3387,3387,3408,3431,3421],label='data',marker='o')
plt.plot(x_y_values_for_912_concise[:,1],label='predicted',marker='o')
plt.xlabel('fiber nr. (counting from the left)')
plt.ylabel('y-position')
plt.legend()
plt.title('for Ar line at 912.58 nm')

x_y_values_for_912_concise[:1]

x_y_values_for_912_concise[:,1]

''

StarcoderdataJupyterScriptsDedupFiltered