# Compute climatologies of SST from OSTIA # Period: Copernicus (1993-2016) # # # Script version 1.0 # # Script works with: # Python version 3.10 # # Package version # numpy: # pandas: # # # # Ver1.0: Created by Mariko Koseki, 03.07.2025 ##--- Import modules ------------------------------## import netCDF4 import numpy as np from numpy import dtype import pandas as pd import math import datetime from datetime import timedelta #from dateutil.relativedelta import relativedelta import os import glob import sys import platform import calendar print('--- Python version ---') print(platform.python_version()) print('--- Package version ---') print('numpy: ', np.__version__) print('pandas: ', pd.__version__) print('-------------------------') print('') ##--- Set path ------------------------------------------## ## Set path to input files ## ''' Set 'ostia_monthly' (path to input NetCDF file) ''' ostia_monthly = '/nird/projects/NS9873K/norcpm/validation/observation/MET/OSTIA/monthly/' ## Set path to output files ## ''' Set 'out_clm_dir' (path to output file) ''' out_dir = '/nird/projects/NS9873K/norcpm/validation/observation/MET/OSTIA/clim/' if not os.path.exists(out_dir): os.makedirs(out_dir) ##--- Set year, month ------------------------------------## ''' Set 'start_year' and 'end_year' Copernicus: 1993-2016 ''' start_year = 1993 # 1993 end_year = 2017 # 2016 y_str_list = ["{}".format(y) for y in range(start_year,end_year)] m_str_list = ["{:02}".format(m) for m in range(1,13)] ''' ##--- Test ------------------------------------## yyyy = '2023' mm = '12' nc_ostia = ostia_monthly + yyyy + '_' + mm + '.nc' print('file name: ', nc_ostia) print('') nc = netCDF4.Dataset(nc_ostia, 'r') ## Check length of variables ## nlon = len(nc.dimensions['lon']) # Longitude: nlat = len(nc.dimensions['lat']) # Latitude: print(nlon, nlat) ## Read variables ## var = nc.variables['sst'][:][:][:] #(time, lat, lon) print(var.shape) #(1, 180, 360) ''' ''' float sst(time, lat, lon) ; sst:standard_name = "sea_surface_foundation_temperature" ; sst:long_name = "analysed sea surface temperature" ; sst:missing_value = -999.f ; sst:comment = " OSTIA foundation SST" ; sst:cell_methods = "time: mean" ; sst:add_offset = 0.f ; sst:scale_factor = 1.f ; sst:units = "Celsius" ; sst:_FillValue = -999.f ; ''' ''' #<> # real_value = (display_value * scale_factor) + add_offset add_offset = nc.variables['sst'].add_offset scale_factor = nc.variables['sst'].scale_factor var_real = (var * scale_factor) + add_offset print(add_offset, scale_factor) print('{:.20g}'.format(var[0,20,20])) print('{:.20g}'.format(var_real[0,20,20])) time = nc.variables['time'][:] #seconds since 1981-01-01 00:00:00 timeunits = nc.variables['time'].getncattr('units') timecalendar = nc.variables['time'].getncattr('calendar') date_sst = netCDF4.num2date(time,timeunits,timecalendar) date_sst = [ f'{i.year}-{i.month:02}' for i in date_sst] print(date_sst) nc.close() sys.exit() ###--- Test end --------------------------------### ''' ## Set length of variables ## nyear = len(y_str_list) nmonth = len(m_str_list) nlat = 180; nlon = 360 months = np.arange(1, nmonth + 1) ##--- Read netCDF files ----------------------------------## ## Calculate climatology ## ostia_clim_3d = np.zeros((nmonth, nlat, nlon)) #(12, 180, 360) for mon_num in range(0, nmonth): mm = m_str_list[mon_num] print('month: ', mm) ##--- Calculate average for nyear ------------------------------------------## ostia_clim_mean = np.zeros((nyear, nlat, nlon)) #(nyear, 180, 360) for year_num in range(0, nyear): yyyy = y_str_list[year_num] print('year: ', yyyy) nc_ostia = ostia_monthly + yyyy + '_' + mm + '.nc' print('file name: ', nc_ostia) print('') nc = netCDF4.Dataset(nc_ostia, 'r') ''' ## Check length of variables ## #nlon = len(nc.dimensions['lon']) # Longitude: #nlat = len(nc.dimensions['lat']) # Latitude: #print(nlon, nlat) ''' ## Read variables ## lons = nc.variables['lon'][:] # Longitude lats = nc.variables['lat'][:] # Latitude var = nc.variables['sst'][:][:][:] #(time, lat, lon)=(1, 180, 360) add_offset = nc.variables['sst'].add_offset scale_factor = nc.variables['sst'].scale_factor var_real = (var * scale_factor) + add_offset ## Convert variables into Numpy ndarray ## var_np = np.array(var_real) ## Add results into a zeros numpy ndarray ## ostia_clim_mean[year_num,:,:] = var_np ## Close netCDF file ## nc.close() mean_2d = np.nanmean(ostia_clim_mean, axis=0) #(180, 360) ostia_clim_3d[mon_num,:,:] = mean_2d print('file shape of 3D', ostia_clim_3d.shape) #(12, 180, 360) print('') print('Calculation completed') print('') print('------------------------') print('') ##--- Save results as a new NetCDF file ------------------------------## print('Save climatology') ## File name ## outfile = out_dir + 'OSTIA_sst_clim_' + str(start_year) + '_' + str(end_year-1) + '.nc' nc2 = netCDF4.Dataset(outfile ,'w', format="NETCDF4") ## Define dimensions ## nc2.createDimension('lat', nlat) nc2.createDimension('lon', nlon) nc2.createDimension('month', nmonth) ## Define variables ## ### latitude ### lat = nc2.createVariable('lat', dtype('float32').char, ('lat',)) lat.units = 'degrees_north' lat.long_name = 'latitude' lat[:] = lats ### longitude ### lon = nc2.createVariable('lon', dtype('float32').char, ('lon',)) lon.units = 'degrees_east' lon.long_name = 'longitude' lon[:] = lons ### month ### month = nc2.createVariable('month', dtype('int32').char, ('month',)) month.units = 'month' month.long_name = 'month' month[:] = months ### climatology ### clim_sst = nc2.createVariable('clim_ostia', dtype('float32').char, ('month', 'lat', 'lon',)) clim_sst.units = 'Celsius' clim_sst.long_name = 'climatology of OSTIA' clim_sst[:,:,:] = ostia_clim_3d ### close NetCDF file ### nc2.close() print('') print('Completed!!')