% initial development : Madlen Kimmritz (madlen.kimmritz@nersc.no), Nov 2019 % This script determines monthly SIV from northern % and southern hemisphere - separately for each member % of reanalysis/hindcast/observation(cont) product and stores the results in a % .mat-file for further use % % model data: siv=hi*parea, as hi= ice volume per unit grid cell area % addpath('/cluster/home/madimm/matlab') %==time for reanalysis yr_start=1950; yr_end =2018; nb_yr=yr_end-yr_start + 1; %================================================================================= %==reanalysis DPATH='/tos-project4/NS9039K/shared/norcpm/cases/NorESM1-CMIP6/'; if (0) %norcpm-cmip6_analysis_19500115, norcpm-cmip6_hindcast (i1 == V1a analysis + hindcasts) CASE='i1'; ASSIM_PATH = [DPATH 'norcpm-cmip6_analysis_19500115/']; HIND_PATH = [DPATH 'norcpm-cmip6_hindcast/norcpm-cmip6_hindcast_']; HINDNAME='norcpm-cmip6_hindcast_' CASENAME='norcpm-cmip6_analysis_19500115_mem'; OUTPUTPATH='/cluster/home/madimm/PROGRAMS_PaperIngoDCPP/PICS/'; else % noresm1-cmip6_analysis_19500115, noresm1-cmip6_hindcast (i2 == V1c analysis + hindcasts) CASE='i2'; ASSIM_PATH = [DPATH 'noresm1-cmip6_analysis_19500115/']; HIND_PATH = [DPATH 'noresm1-cmip6_hindcast/noresm1-cmip6_hindcast_']; HINDNAME='noresm1-cmip6_hindcast_' OUTPUTPATH='/cluster/home/madimm/PROGRAMS_PaperIngoDCPP/PICS/'; CASENAME='noresm1-cmip6_analysis_19500115_mem'; end; SIV_Arctic=zeros(30,nb_yr*12); SIV_SOcean=zeros(30,nb_yr*12); %=== GRID idm=320; jdm=384; grid_path='/cluster/home/madimm/norcpm_link/SFE/'; plon=ncgetvar([grid_path 'grid.nc'],'plon'); plat=ncgetvar([grid_path 'grid.nc'],'plat'); parea=ncgetvar([grid_path 'grid.nc'],'parea'); areaN=parea*1.;areaN(:,1:180)=0.; %use upper half areaS=parea*1.;areaS(:,180:end)=0.; %======= SIV (North and SOUTH) for reanalysis if (1) for mem=1:1:30 mem time=0; for year=yr_start:1:yr_end year for mm=1:12 time=time+1; filename=[ASSIM_PATH CASENAME num2str(mem,'%2.2d') '/ice/hist/' CASENAME num2str(mem,'%2.2d') '.cice.h.' num2str(year) '-' num2str(mm,'%2.2d') '.nc']; tmp=ncgetvar(filename, 'hi'); %unit: m, sea ice volume per unit area AN=nansum(nansum(areaN.*tmp)); AS=nansum(nansum(areaS.*tmp)); SIV_Arctic(mem,time)=AN(1,1); SIV_SOcean(mem,time)=AS(1,1); end; end; end; save([OUTPUTPATH '/SIV_' CASE '_time.mat'], 'SIV_Arctic', 'SIV_SOcean'); end; %================================================================================= %==== FREE: /projects/NS9039K/shared/norcpm/cases/NorESM1-CMIP6/noresm1-cmip6_historical (free historical ensemble) FREE_PATH = [DPATH 'noresm1-cmip6_historical_19500101/']; FREE_PATH1 = [DPATH 'noresm1-cmip6_historical_20150101/']; OUTPUTPATH='/cluster/home/madimm/PROGRAMS_PaperIngoDCPP/PICS/'; CASENAMEF='noresm1-cmip6_historical_19500101_mem'; CASENAMEF1='noresm1-cmip6_historical_20150101_mem'; nbf=2014-1950+1; SIV_Arctic_FREE=zeros(30,nb_yr*12); SIV_SOcean_FREE=zeros(30,nb_yr*12); for mem=1:1:30 mem time=0; for year=yr_start:1:2014 year for mm=1:12 time=time+1; filename=[FREE_PATH CASENAMEF num2str(mem,'%2.2d') '/ice/hist/' CASENAMEF num2str(mem,'%2.2d') '.cice.h.' num2str(year) '-' num2str(mm,'%2.2d') '.nc']; tmp=ncgetvar(filename, 'hi'); %unit: m, sea ice volume per unit area AN=nansum(nansum(areaN.*tmp)); AS=nansum(nansum(areaS.*tmp)); SIV_Arctic_FREE(mem,time)=AN(1,1); SIV_SOcean_FREE(mem,time)=AS(1,1); end; end; %== for year=2015:1:yr_end year for mm=1:12 time=time+1; filename=[FREE_PATH1 CASENAMEF1 num2str(mem,'%2.2d') '/ice/hist/' CASENAMEF1 num2str(mem,'%2.2d') '.cice.h.' num2str(year) '-' num2str(mm,'%2.2d') '.nc']; tmp=ncgetvar(filename, 'hi'); %unit: m, sea ice volume per unit area AN=nansum(nansum(areaN.*tmp)); AS=nansum(nansum(areaS.*tmp)); SIV_Arctic_FREE(mem,time)=AN(1,1); SIV_SOcean_FREE(mem,time)=AS(1,1); end; end; end; save([OUTPUTPATH '/SIV_FREE_time.mat'], 'SIV_Arctic_FREE', 'SIV_SOcean_FREE'); %================================================================================= %======= SIA (North and SOUTH) for observations if(0) %=== PIOMAS DATA for ARCTIC %PIOMAS MNTHLY DATA (unit: 1000 km^3) SIV_PIOMAS_Arctic=zeros(2018-1979+1,13); SIV_PIOMAS_Arctic(1,1:12) = [27.704 30.171 32.045 32.951 32.295 29.787 23.661 18.409 16.911 17.854 20.124 23.201]; %1979 SIV_PIOMAS_Arctic(2,1:12) = [26.530 29.154 31.128 32.235 31.807 29.152 22.865 17.790 16.316 17.325 19.424 22.451]; SIV_PIOMAS_Arctic(3,1:12) = [25.317 27.770 29.818 30.744 30.034 26.817 20.173 14.635 12.808 13.965 16.145 19.115]; SIV_PIOMAS_Arctic(4,1:12) = [ 22.633 25.486 27.653 28.968 28.302 25.616 19.665 14.834 13.508 14.911 17.766 21.066]; SIV_PIOMAS_Arctic(5,1:12) = [ 24.457 27.279 29.385 30.386 30.170 27.758 21.854 16.738 15.200 16.427 18.954 21.962]; SIV_PIOMAS_Arctic(6,1:12) = [ 24.733 27.209 29.151 30.330 29.816 27.093 20.941 16.148 14.635 15.664 18.097 21.292]; SIV_PIOMAS_Arctic(7,1:12) = [ 24.639 27.264 29.453 30.864 30.517 27.488 20.768 16.012 14.583 15.838 18.289 21.285]; SIV_PIOMAS_Arctic(8,1:12) = [ 24.678 27.518 29.770 30.934 30.456 27.861 21.953 17.293 16.076 17.378 19.857 22.738]; SIV_PIOMAS_Arctic(9,1:12) = [ 26.025 28.840 30.663 31.790 31.502 28.814 22.274 16.708 15.353 16.703 19.255 22.299]; SIV_PIOMAS_Arctic(10,1:12) = [ 25.471 28.132 30.283 31.200 30.213 27.230 21.098 16.328 14.987 16.246 18.925 22.125]; SIV_PIOMAS_Arctic(11,1:12) = [ 25.363 27.837 29.424 30.112 29.621 27.122 21.075 16.198 14.770 15.941 18.654 21.797]; SIV_PIOMAS_Arctic(12,1:12) = [ 24.967 27.632 29.393 29.907 28.905 25.494 19.515 14.885 13.815 15.218 18.258 21.522]; SIV_PIOMAS_Arctic(13,1:12) = [ 24.689 27.477 29.665 30.742 30.005 26.792 20.014 15.005 13.593 15.017 17.688 20.894]; SIV_PIOMAS_Arctic(14,1:12) = [ 24.183 26.906 28.680 29.642 29.259 26.849 20.899 15.915 15.078 16.580 19.157 22.328]; SIV_PIOMAS_Arctic(15,1:12) = [ 25.300 27.720 29.462 30.427 29.575 25.805 18.509 13.479 12.441 13.964 16.983 20.398]; SIV_PIOMAS_Arctic(16,1:12) = [ 23.853 26.673 28.714 29.735 29.312 26.552 20.030 14.895 13.855 15.196 17.940 21.211]; SIV_PIOMAS_Arctic(17,1:12) = [ 24.199 26.431 27.953 28.447 27.311 23.810 17.186 12.438 11.236 12.128 14.837 18.178]; SIV_PIOMAS_Arctic(18,1:12) = [ 21.648 24.474 26.387 27.445 27.255 24.866 19.282 14.896 13.949 15.394 17.398 20.152]; SIV_PIOMAS_Arctic(19,1:12) = [ 23.497 26.315 28.317 29.370 28.741 25.716 19.242 14.381 13.228 14.199 16.615 19.916]; SIV_PIOMAS_Arctic(20,1:12) = [ 23.472 26.355 28.300 29.409 28.823 25.361 18.577 13.286 11.626 12.812 15.640 18.951]; SIV_PIOMAS_Arctic(21,1:12) = [ 22.428 25.324 27.299 28.450 27.924 24.881 18.517 13.076 11.045 12.433 15.377 18.394]; SIV_PIOMAS_Arctic(22,1:12) = [ 21.712 24.297 26.174 27.153 26.689 23.881 17.429 12.422 11.084 12.313 14.901 17.952]; SIV_PIOMAS_Arctic(23,1:12) = [ 21.259 24.128 26.414 27.632 26.870 23.828 17.901 13.397 12.275 13.321 15.764 18.605]; SIV_PIOMAS_Arctic(24,1:12) = [ 21.837 24.793 26.655 27.428 26.822 23.685 17.164 12.127 10.846 11.890 14.587 17.844]; SIV_PIOMAS_Arctic(25,1:12) = [ 21.236 24.127 26.308 27.244 26.273 22.911 16.486 11.496 10.284 11.204 13.686 16.848]; SIV_PIOMAS_Arctic(26,1:12) = [ 20.050 22.655 24.809 25.750 25.308 22.706 16.577 11.515 10.037 11.247 13.966 17.238]; SIV_PIOMAS_Arctic(27,1:12) = [ 20.254 22.667 24.820 26.045 25.354 21.617 15.225 10.708 9.283 10.180 12.847 15.983]; SIV_PIOMAS_Arctic(28,1:12) = [ 19.307 22.036 24.051 25.105 24.337 20.904 14.692 10.404 9.113 9.810 12.263 15.039]; SIV_PIOMAS_Arctic(29,1:12) = [ 18.322 20.813 22.995 23.752 23.094 19.190 12.119 7.615 6.526 7.107 10.404 14.167]; SIV_PIOMAS_Arctic(30,1:12) = [ 18.515 21.527 23.823 24.983 24.122 20.576 14.154 9.201 7.250 8.267 11.700 15.134]; SIV_PIOMAS_Arctic(31,1:12) = [ 18.794 21.682 23.794 24.944 23.880 19.732 12.833 8.272 6.925 7.630 10.759 14.184]; SIV_PIOMAS_Arctic(32,1:12) = [ 17.667 20.581 23.082 24.097 22.223 17.142 10.244 5.927 4.742 6.198 9.480 12.925]; SIV_PIOMAS_Arctic(33,1:12) = [ 16.206 19.320 21.390 22.500 21.139 16.500 9.546 5.489 4.480 5.719 9.250 12.989]; SIV_PIOMAS_Arctic(34,1:12) = [ 16.886 19.585 21.924 23.119 21.714 16.002 9.264 4.964 3.787 5.001 8.217 12.133]; SIV_PIOMAS_Arctic(35,1:12) = [ 15.799 19.318 21.964 23.122 21.865 17.538 10.538 6.394 5.479 6.953 10.076 13.787]; SIV_PIOMAS_Arctic(36,1:12) = [ 17.406 19.850 21.802 22.935 21.913 17.680 11.951 8.150 6.972 8.159 11.481 15.074]; SIV_PIOMAS_Arctic(37,1:12) = [ 18.445 21.462 23.213 24.230 23.027 18.557 11.649 7.086 5.852 7.013 10.298 14.001]; SIV_PIOMAS_Arctic(38,1:12) = [ 17.185 19.592 21.524 22.459 21.026 16.493 10.257 5.941 4.530 5.511 7.834 11.206]; SIV_PIOMAS_Arctic(39,1:12) = [ 14.641 17.366 19.576 20.664 19.834 15.400 9.141 5.450 4.679 6.108 9.179 12.638]; SIV_PIOMAS_Arctic(40,1:12) = [ 16.022 18.639 20.844 22.245 21.429 17.206 10.295 6.203 5.081 5.776 9.426 13.172]; %2018 for i=1:40 SIV_PIOMAS_Arctic(i,13) =mean(SIV_PIOMAS_Arctic(i,1:12)); end; SIV_PIOMAS_Arctic= SIV_PIOMAS_Arctic(:,13); save([OUTPUTPATH '/SIV_PIOMAS_time.mat'], 'SIV_PIOMAS_Arctic'); %========================================================== %=== combined data from cryosat and hadisst1 %==OBS + OBS grid det_obs = 1 if (det_obs) yr_starto=2002; yr_endo=2017; nobs=2017-2002+1; %=========== SIT DATA %SIT availability: (2002-2017), cryosat2 and envisat OSIT_PATH='/cluster/shared/noresm/norcpm/Obs/ICET/Merged/'; ipivSIT=int16(ncread('/cluster/shared/noresm/norcpm/Input/NorESM/f19_g16/ICET/pivots_ICET.nc','ipiv')); jpivSIT=int16(ncread('/cluster/shared/noresm/norcpm/Input/NorESM/f19_g16/ICET/pivots_ICET.nc','jpiv')); dtlon=432; dtlat=432; tmpto=zeros(dtlon,dtlat); SIT_Arctic=zeros(dtlon,dtlat,nobs*12); SIC_Arctic=zeros(dtlon,dtlat,nobs*12); SIT_found=zeros(1,nobs*12); %1 when data available, 0 else, availability in winter months, we take average only over %those months in year, where data available (= good estimate of annual avg) time=0; for year=yr_starto:1:yr_endo year for mm=1:12 time=time+1; if exist([OSIT_PATH num2str(year) '_' num2str(mm,'%2.2d') '.nc'], 'file') == 2 SIT_found(1,time)=1; time tmpto=ncgetvar([OSIT_PATH num2str(year) '_' num2str(mm,'%2.2d') '.nc'], 'sea_ice_thickness'); tmpco=ncgetvar([OSIT_PATH num2str(year) '_' num2str(mm,'%2.2d') '.nc'], 'sea_ice_concentration'); tmpco=tmpco./100.; SIC_Arctic(:,:,time)= tmpco; SIT_Arctic(:,:,time)= tmpto; % end; end; end; %========= SIC DATA %SIC availability: (1870-2018), HadISST1 OSIC_PATH='/cluster/projects/nn9039k/NorCPM/Obs/ICEC/HADISST2/HadISST_ice_1870_2018.nc'; ipivSIC=int16(ncread('/cluster/projects/nn9039k/NorCPM/Old/Input/NorESM/f19_g16/HADISST2/pivots_hadi-micom.nc','ipiv')); jpivSIC=int16(ncread('/cluster/projects/nn9039k/NorCPM/Old/Input/NorESM/f19_g16/HADISST2/pivots_hadi-micom.nc','jpiv')); SIC_HadISST_Arctic=zeros(idm,jdm,nobs*12); SIC_HadISST_Arctic_T=zeros(dtlon,dtlat,nobs*12); area_T=zeros(dtlon,dtlat); dlon=360; dlat=180; tmpo=zeros(dlon,dlat); %1870 - 2001: 132 years (132*12)+1 = Jan 2002 time1=(132*12)+1; obsice = ncgetvar(OSIC_PATH,'sic'); obsice=obsice(:,:,time1:end); d=size(obsice); time=0; for year=yr_starto:1:yr_endo year for mm=1:12 time=time+1; tmpo=obsice(:,:,time); %== transfer to model grid nb =zeros(idm,jdm); for i=1:idm for j=1:jdm SIC_HadISST_Arctic(i,j,time)=SIC_HadISST_Arctic(i,j,time) + tmpo(ipivSIC(i,j),jpivSIC(i,j)); nb(i,j)=nb(i,j)+1; end end SIC_HadISST_Arctic(:,:,time)= SIC_HadISST_Arctic(:,:,time)./nb; %==== nbt =zeros(dtlon,dtlat); for i=1:dtlon for j=1:dtlat SIC_HadISST_Arctic_T(i,j,time)=SIC_HadISST_Arctic_T(i,j,time) + SIC_HadISST_Arctic(ipivSIT(i,j),jpivSIT(i,j),time); nbt(i,j)=nbt(i,j)+1; end end SIC_HadISST_Arctic_T(:,:,time)= SIC_HadISST_Arctic_T(:,:,time)./nbt; end; end; %=== nbt =zeros(dtlon,dtlat); for i=1:dtlon for j=1:dtlat area_T(i,j)=area_T(i,j) + areaN(ipivSIT(i,j),jpivSIT(i,j)); nbt(i,j)=nbt(i,j)+1; end end area_T(:,:)= area_T(:,:)./nbt; %======== SIV DATA time=0; SIV_obsmerged_Arctic=zeros(dtlon,dtlat,nobs*12); for year=yr_starto:1:yr_endo year for mm=1:12 time=time+1; % SIV_obsmerged_Arctic(:,:,time)=SIC_HadISST_Arctic_T(:,:,time).*SIT_Arctic(:,:,time); SIV_obsmerged_Arctic(:,:,time)=SIC_Arctic(:,:,time).*SIT_Arctic(:,:,time); % SIV_obsmerged_Arctic(:,:,time)=(SIT_Arctic(:,:,time)); SIV_obsmerged_Arctic(:,:,time)=SIV_obsmerged_Arctic(:,:,time).*area_T; end; end; SIV_obsmerged_Arctic_time=zeros(1,nobs*12); time=0; for year=yr_starto:1:yr_endo year for mm=1:12 time=time+1; SIV_obsmerged_Arctic_time(1,time)=nansum(nansum(SIV_obsmerged_Arctic(:,:,time))); end; end; SIV_obsmerged_Arctic_time_merged = SIV_obsmerged_Arctic_time.*1.; SIV_obsmerged_Arctic_time_HadISST = SIV_obsmerged_Arctic_time.*1.; save([OUTPUTPATH '/SIV_OBSmerged_time.mat'], 'SIV_obsmerged_Arctic_time_merged', 'SIV_obsmerged_Arctic_time_HadISST', 'SIT_found'); end; end; %================================================================================= % == hindcast (10 members, 10 years hindcast length) yr_start=1960; yr_end =2018; nb_yr=yr_end-yr_start + 1; SIV_Arctic_hind=zeros(10,nb_yr, 10*12); SIV_SOcean_hind=zeros(10,nb_yr, 10*12); for mem=1:10 nyear=0; for year=yr_start:1:yr_end nyear=nyear+1; CASENAME=[HINDNAME num2str(year) '1015_mem' num2str(mem,'%2.2d')] for mm=1:120 yr1=year + floor(((mm+10-1))./12); %lead month and year mm1=mod((mm+10-1),12)+1; %disp([year, mm, mm1, yr1]); filename=[HIND_PATH num2str(year) '1015/' CASENAME '/ice/hist/' CASENAME '.cice.h.' num2str(yr1) '-' num2str(mm1,'%2.2d') '.nc']; %disp(filename); tmp=ncgetvar(filename, 'hi'); AN=0.;AS=0.; AN=nansum(nansum(areaN.*tmp)); AS=nansum(nansum(areaS.*tmp)); SIV_Arctic_hind(mem,nyear,mm)=AN(1,1); SIV_SOcean_hind(mem,nyear,mm)=AS(1,1); end; end; end; save([OUTPUTPATH '/SIV_' CASE 'hind_time.mat'], 'SIV_Arctic_hind', 'SIV_SOcean_hind');