# Compare sea ice area/extent in the Arctic # 1. OISST (NOAA; ver2, ver2.1) # 2. OSTIA (MET/Coperinicus) # # # Plot sea ice area/extent # Compare observational data # # # Script version 1.0 # # Python version: 3.10 # # # Ver1.0: Created by Mariko Koseki, 18.06.2025 import numpy as np import pandas as pd import datetime from datetime import timedelta from dateutil.relativedelta import relativedelta import statistics import math import collections import matplotlib import matplotlib.pyplot as plt import matplotlib.dates as mdates import matplotlib.ticker as ptick 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('matplotlib: ', matplotlib.__version__) ##--- Input -----------------------------------------## args = sys.argv if len(args) == 2: year_start = args[1] if len(year_start) == 4: yyyy1 = int(year_start[0:4]) yyyy2 = yyyy1 else: print('---How to use this script---') print('python', args[0], '') print('or') print('python', args[0], ' ') print(' = start year; = end year') print('') print('Example: python', args[0], '2023') print('Example: python', args[0], '1995 2023') print('') sys.exit() elif len(args) == 3: year_start = args[1] year_end = args[2] if (len(year_start) == 4) & (len(year_end) == 4): yyyy1 = int(year_start[0:4]) yyyy2 = int(year_end[0:4]) else: print('') print('---How to use this script---') print('python', args[0], '') print('or') print('python', args[0], ' ') print(' = start year; = end year') print('') print('Example: python', args[0], '2023') print('Example: python', args[0], '1995 2023') print('') sys.exit() else: print('') print('---How to use this script---') print('python', args[0], '') print('or') print('python', args[0], ' ') print(' = start year; = end year') print('') print('Example: python', args[0], '2023') print('Example: python', args[0], '1995 2023') print('') sys.exit() ##--- Create list of year ------------------------------------------------------## yyyy_int_list = [] for yyyy in range(yyyy1, yyyy2+1): yyyy_int_list.append(yyyy) print(yyyy_int_list) ##--- Set path ------------------------------------------## ## Input ## ''' Set path to input directory ### ''' seaice_dir = '/nird/projects/NS9873K/norcpm/validation/seaicepre/' ### OISST ### #out_ice_oisst21 = seaice_dir + 'oisst/' + ver21 + '/' #out_ice_oisst2 = seaice_dir + 'oisst/' + ver2 + '/' out_oisst_post = seaice_dir + 'oisst/merge/' ### OSTIA ### out_ostia = seaice_dir + 'ostia/yearly/' ## Output ## ''' Set path to output file ''' out_dir = '/nird/projects/NS9873K/www/users/mariko/fig/seaice/compare/' if not os.path.exists(out_dir): os.makedirs(out_dir) for YEAR in yyyy_int_list: print('') print('Year: ', YEAR) print('') ##--- Read CSV files ------------------------------## ## Paths to CSV files ## file_oisst = out_oisst_post + 'seaice_oisst_merge_' + str(YEAR) + '.csv' file_ostia = out_ostia + 'seaice_ostia_merge_' + str(YEAR) + '.csv' print('') print('file of OISST: ', file_oisst) print('') print('file of OSTIA: ', file_ostia) print('') ## CSV files -> pandas.DataFrame ## ### OISST ### ice_oisst = pd.read_csv(file_oisst) #print(ice_oisst.tail(50)) ### OSTIA ### ice_ostia = pd.read_csv(file_ostia) #print(ice_ostia) ## Change 'object' to 'datetime' ## ### OISST ### ice_oisst['dtime'] = pd.to_datetime(ice_oisst['dtime']) ### OSTIA ### ice_ostia['dtime'] = pd.to_datetime(ice_ostia['dtime']) ##--- Plot interpolated sea ice area ------------------------## df_oisst = ice_oisst.copy() df_ostia = ice_ostia.copy() #date = df_oisst.iloc[0,0] #yyyy = date.year time_obs = df_oisst['dtime'] #print(len(time_obs)) ### AREA ### area_oisst = df_oisst['sea_ice_area'] area_inter_oisst = df_oisst['sea_ice_area_interpolated'] area_ostia = df_ostia['sea_ice_area'] area_inter_ostia = df_ostia['sea_ice_area_interpolated'] ### EXTENT ### extent_oisst = df_oisst['sea_ice_extent'] extent_inter_oisst = df_oisst['sea_ice_extent_interpolated'] extent_ostia = df_ostia['sea_ice_extent'] extent_inter_ostia = df_ostia['sea_ice_extent_interpolated'] fig3 = plt.figure(figsize=(8,5)) ax3 = fig3.add_subplot(111) ### AREA ### ax3.plot(time_obs, area_inter_oisst, linewidth = 4.7, label= 'Sea ice area (OISST/interpolated/ver2.1&ver2)', alpha=0.7, color='b') ax3.plot(time_obs, area_oisst, linewidth = 1.3, label = 'Sea ice area(OISST/ver2.1)', color='k') ax3.plot(time_obs, area_inter_ostia, linewidth = 4.7, label = 'Sea ice area(OSTIA/interpolated)', alpha=0.7, color='magenta') ax3.plot(time_obs, area_ostia, linewidth = 1.3, label = 'Sea ice area(OSTIA)', color='r') ### EXTENT ### #ax3.plot(time_obs, extent_inter_oisst, linewidth = 4.7, label= 'Sea ice extent (interpolated)', color='magenta', alpha=0.7) #ax3.plot(time_obs, extent_oisst, linewidth = 1.3, label = 'Sea ice extent', color='r') ## Set xlim, ylim ## span = [datetime.datetime(YEAR, 1, 1), datetime.datetime(YEAR, 12, 31)] ax3.set_xlim(span) ax3.set_ylim([0, 16000000]) ### Set X axis ### ax3.set_xlabel('') ax3.xaxis.set_major_locator(mdates.MonthLocator(interval=1)) # interval: 1 month plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b')) ### Set Y axis ### ax3.yaxis.set_major_formatter(ptick.ScalarFormatter(useMathText=True)) ax3.ticklabel_format(style="sci", axis="y", scilimits=(6,6)) # 1*10^6 ## Add grid ## ax3.minorticks_on() ax3.grid(color = "gray", linestyle = "--") ## Add legend ## ax3.legend(fontsize = 'small') ## Add labels, title ## ax3.set_ylabel('Sea Ice Area [$\mathrm{km^2}$]', fontsize = 12) fig_title = 'Daily Arctic sea ice area (' + str(YEAR) + ')' plt.suptitle(fig_title, fontsize = 15) ## Save figures ## fig_name = out_dir + '/sea_ice_area_interpolated_' + str(YEAR) + '.png' fig3.savefig(fig_name, dpi=200, format='png') plt.show() print('') print('Completed!')