<jupyter_start><jupyter_text>This makes images from the goes collocations<jupyter_code>%matplotlib inline
from matplotlib import style
style.use('ggplot')
import numpy as np
import matplotlib.pyplot as plt
import numpy.ma as ma
from netCDF4 import Dataset, date2index, num2date
from palettable import colorbrewer
from datetime import datetime, timedelta
import datetime as dt
from math import atan2
import xarray as xr

#make maps of daily goes ssts
import matplotlib as mpl
#import openpyxl
from mpl_toolkits.basemap import Basemap
#from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
from copy import copy
import matplotlib.pyplot as plt
import matplotlib.colors as colors
import matplotlib.mlab as mlaba
palette = copy(plt.cm.jet)
palette.set_over('r', 1.0)
palette.set_under('g', 1.0)
palette.set_bad(alpha = 0.0)

filename_goes_sst = 'F:/data/cruise_data/saildrone/baja-2018/goes_sst_collocation_data.npy'
filename_usv='F:/data/cruise_data/saildrone/baja-2018/saildrone-gen_4-baja_2018-EP-sd1002-ALL-1_min-v1.nc'
#filename_goes_sst = './goes_sst_collocation_data.npy'
#filename_usv='./saildrone-gen_4-baja_2018-EP-sd1002-ALL-1_min-v1.nc'

dataset = xr.open_dataset(filename_usv)
ilen=(len(dataset.LAT['obs']))
print(dataset.LAT[0,1].values*1)
print(type(dataset))
lats_usv=dataset.LAT[0,:].values
lons_usv=dataset.LON[0,:].values
date_1970 = dt.datetime(1970,1,1,0,0,0) # start date is1/1/1970
skin=dataset.TEMP_IR_MEAN[0,:].values
tair=dataset.TEMP_AIR_MEAN[0,:].values
bulk1=dataset.TEMP_CTD_MEAN[0,:].values
bulk2=dataset.TEMP_O2_MEAN[0,:].values
uspd=dataset.UWND_MEAN[0,:].values
vspd=dataset.VWND_MEAN[0,:].values
wing=dataset.HDG_WING[0,:].values
wing_ang=dataset.WING_ANGLE[0,:].values
pitch=dataset.PITCH[0,:].values
heading=dataset.HDG_WING[0,:].values
yaw_heading=dataset.HDG[0,:].values
roll=dataset.ROLL[0,:].values
#sec_since_1970 = dataset.TIME[0,:].values #'seconds since 1970-01-01T00:00:00+00:00'
tdim=len(skin)
print('here')
tem_date=[0]*tdim #np.zeros(tdim)
tem_dy_from=np.zeros(tdim)
tem_lmt=np.zeros(tdim)
#tem_date=[0]*tdim
#tem_dy_from=[0]*tdim
#tem_lmt=[0]*tdim
wdir=np.zeros(tdim)
print(dataset.TIME.shape)
#for i in range(0,10):
#    print(float(dataset.TIME[0,i].values)/86400000000000.)
for i in range(0,tdim):
    tem_dy=float(dataset.TIME[0,i].values)/86400000000000.
    tem_dy_from[i]=float(dataset.TIME[0,i].values)/86400000000000.-float(dataset.TIME[0,0].values)/86400000000000.
    tem_date[i]=date_1970+dt.timedelta(days=tem_dy)  #create new time array that can be queried for year etc
    tem_lmt[i]=tem_date[i].hour+(dataset.LON[0,i].values/360.)*24.
    if tem_lmt[i]>24:
        tem_lmt[i]=tem_lmt[i]-24.
    if tem_lmt[i]<0:
        tem_lmt[i]=tem_lmt[i]+24.
    wdir[i]=atan2(uspd[i],vspd[i])*180/3.14159
    if i<5:
        print(tem_date[i],wdir[i])
        print(tem_dy_from[i])
wspd=(uspd**2+vspd**2)**.5
print('time:', tem_dy_from[1:5],tem_dy_from[-5:-1])

tem_dyr=np.zeros(tdim)
for i in range(0,tdim):
    tem_dyr[i] = int(tem_date[i].timetuple().tm_yday)


print(tdim)
adir_goes='F:/data/cruise_data/saildrone/baja-2018/figs/goes_sst_daily/'    
#adir_goes='./figs/goes_sst_daily/'    
gsst = np.load(filename_goes_sst)    

tdim=gsst.shape[0]
istart=0
lat0=59.974998474121094
lon0=134.97500610351562
dx=0.05
dy=0.05
gsst=np.zeros(tdim)
#gsst_box=np.zeros(tdim,10,10)
#print(0*dx-lon0,1*dx-lon0,2*dx-lon0,2399*dx-lon0)
#print(0*dy-lat0,1*dy-lat0,2*dy-lat0,2399*dy-lat0)
inew_data=1
save_date=datetime(2000,1,1)
iy=int(round((lats_usv.mean()+lat0)/dy))
ix=int(round((lons_usv.mean()+lon0)/dx))   

for i in range(0,tdim):
    if tem_date[i].minute!=0:
        continue
#if new day set new center point
    if tem_date[i].day!=save_date.day:
        print('here')
        save_date=tem_date[i]
#make daily map
        if inew_data==0:
            print('here2')
            mask = (sst_new < -10)
            sst_new[mask]=-9999.
            sst_new2 = sst_new[0,:,:]
            sstx = ma.masked_values (sst_new2, -9999.)
            sst5x=ma.swapaxes(sstx,0,1)
            cmap = colorbrewer.get_map('Spectral', 'diverging', 11, reverse=True).mpl_colormap
            sub_lons=lons_usv[tem_dyr==day_of_year]
            sub_lats=lats_usv[tem_dyr==day_of_year]
            min_lon=min(sub_lons)-5
            max_lon=min(sub_lons)+5
            min_lat=min(sub_lats)-5
            max_lat=min(sub_lats)+5
            center_lon=.5*(max_lon-min_lon)+min_lon
            center_lat=.5*(max_lat-min_lat)+min_lat
            print(min_lon,max_lon,min_lat,max_lat,center_lon,center_lat)

            fig = plt.figure(figsize=(11.7,8.3))
            m = Basemap(projection='merc', lat_0 = center_lat, lon_0 = center_lon,    resolution = 'l', area_thresh = 0.1,    
                        llcrnrlon=min_lon, llcrnrlat=min_lat,urcrnrlon=max_lon, urcrnrlat=max_lat)
            m.bluemarble()
            lat_grid, lon_grid = np.meshgrid(lat, lon)
            x,y = m(lon_grid,lat_grid)
            cs = m.pcolormesh(x,y, sst5x - 273.15, cmap=cmap, vmin = 12, vmax = 19)
            x,y = m(sub_lons,sub_lats)
            m.plot(x,y,'w', linewidth=3)
            m.plot(x,y,color='magenta', linewidth=1)
            m.colorbar(location="bottom",label='SST [C]') # draw colorbar
            plt.show()
            filename_png=adir_goes+'goes_cruise_sst_'+str(day_of_year)+'_zoom.png'
            fig.savefig(filename_png, transparent=False, format='png')
            print('done') 
            mask = sst_new < 500
            sst_new[mask] = -9999
            mask = (sst_new < -10)
            sst_new[mask]=-9999.
            sst_new2 = sst_new[0,:,:]
            sstx = ma.masked_values (sst_new2, -9999.)
            sst5x=ma.swapaxes(sstx,0,1)
            cmap = colorbrewer.get_map('Spectral', 'diverging', 11, reverse=True).mpl_colormap
            sub_lons=lons_usv[tem_dyr==day_of_year]
            sub_lats=lats_usv[tem_dyr==day_of_year]
            min_lon=min(sub_lons)-5
            max_lon=min(sub_lons)+5
            min_lat=min(sub_lats)-5
            max_lat=min(sub_lats)+5
            center_lon=.5*(max_lon-min_lon)+min_lon
            center_lat=.5*(max_lat-min_lat)+min_lat
            print(min_lon,max_lon,min_lat,max_lat,center_lon,center_lat)

            fig = plt.figure(figsize=(11.7,8.3))
            m = Basemap(projection='merc', lat_0 = center_lat, lon_0 = center_lon,    resolution = 'l', area_thresh = 0.1,    
                        llcrnrlon=min_lon, llcrnrlat=min_lat,urcrnrlon=max_lon, urcrnrlat=max_lat)
            m.bluemarble()
            lat_grid, lon_grid = np.meshgrid(lat, lon)
            x,y = m(lon_grid,lat_grid)
            cs = m.pcolormesh(x,y, sst5x - 273.15, cmap=cmap, vmin = 12, vmax = 19)
            x,y = m(sub_lons,sub_lats)
            m.plot(x,y,'w', linewidth=3)
            m.plot(x,y,color='magenta', linewidth=1)
            m.colorbar(location="bottom",label='SST [C]') # draw colorbar
            plt.show()
   
    syr=str(tem_date[i].year).zfill(4)
    smon=str(tem_date[i].month).zfill(2)
    sdym=str(tem_date[i].day).zfill(2)
    shr=str(tem_date[i].hour).zfill(2)
    day_of_year=tem_date[i].timetuple().tm_yday
    sday_of_year=str(tem_date[i].timetuple().tm_yday).zfill(3)
    if tem_date[i].hour==0:
        d=tem_date[i]+ timedelta(days=-1) #increment day by -1 if 0 hr
        sday_of_year=str(d.timetuple().tm_yday).zfill(3)
    adir_goes='https://opendap.jpl.nasa.gov/opendap/OceanTemperature/ghrsst/data/GDS2/L3C/AMERICAS/GOES16/OSISAF/v1/'
    fgoes='0000-OSISAF-L3C_GHRSST-SSTsubskin-GOES16-ssteqc_goes16_'
    dstr=syr+smon+sdym+shr
    dstr2=syr+smon+sdym+'_'+shr
    file_goes=adir_goes+syr+'/'+sday_of_year+'/'+ dstr + fgoes +dstr2+'0000-v02.0-fv01.0.nc'   
    r = requests.get(file_goes)
    if r.status_code != requests.codes.ok:
        continue
    print(r.status_code,file_goes)
    nc = Dataset(file_goes)
    sst = nc.variables['sea_surface_temperature'][:,iy-200:iy+200,ix-200:ix+200]
    nc.close()
#if first time, initialize new saved masked array read in lat/lons
    if inew_data==1:
        inew_data = 0
        nc = Dataset(file_goes)
        sst = nc.variables['sea_surface_temperature'][:,iy-200:iy+200,ix-200:ix+200]
        sst_new = ma.copy(sst)
        lat = nc.variables['lat'][iy-200:iy+200]
        lon = nc.variables['lon'][ix-200:ix+200]
        nc.close()
        save_date=tem_date[i]
        print('init data',inew_data)
#add new data to saved masked array
    #print(sst.shape,sst_new.shape)
    mask = sst_new < -10
    sst_new[mask] = sst[mask]
    

    

adir_goes='F:/data/cruise_data/saildrone/baja-2018/figs/goes_collocations/'    
#adir_goes='./figs/goes_collocations/'    

print(len(tem_dy_from))
print(len(skin-bulk1),tdim)
print(len(gsst))
item=np.zeros(tdim)
for i in range(0,tdim):
    item[i]=int(round(tem_dy_from[i]/5.))
print(item[0:5000:500],max(item),item[tdim-1])
print(tem_dy_from[0:5000:500],tem_dy_from[tdim-1])
#dt=gsst-bulk1
i=0
mindy=min(tem_date)
maxdy=max(tem_date)
print(maxdy-mindy)
intdy=(maxdy-mindy)/12.
fig, ax1 = plt.subplots() 
plt.scatter(tem_date,skin-bulk1,c='b')
plt.scatter(tem_date,gsst-273.15-bulk1,c='r')
plt.xlabel('Date')
plt.ylabel('GOES SST - T$_{ctd}$')
for i in range(0,12):
    plt.ylim(-3, 3)
    plt.xlim(mindy+i*intdy,mindy+(i+1)*intdy)
    #plt.show()
    fig_fname=adir_goes+'skinbulk_goes_TS' + str(i) +'.png'
    print(fig_fname)
    fig.savefig(fig_fname, transparent=False, format='png')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.15543822597676873, 'fraction_numerical': 0.05385427666314678, 'mean_word_length': 3.6426470588235293, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 6, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 4, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '1', 'dir': 'jupyter-structured-clean-dedup', 'id': '254122', 'license': 'Apache-2.0', 'n_tokens_mistral': 4351, 'n_tokens_neox': 4075, 'n_words': 501, 'path': '.ipynb_checkpoints/goes16_collocation-images-checkpoint.ipynb', 'repo_name': 'riddhimap/Saildrone', 'size': '9470'}