From c2bde80f7689dd5e329e08d93e757372e3f7334d Mon Sep 17 00:00:00 2001 From: turner Date: Mon, 20 Nov 2017 20:21:17 +0000 Subject: [PATCH 1/8] Add new script for the CICE t-test to test for non-BFB --- t-test.py | 221 ++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 221 insertions(+) create mode 100755 t-test.py diff --git a/t-test.py b/t-test.py new file mode 100755 index 000000000..6230f76d5 --- /dev/null +++ b/t-test.py @@ -0,0 +1,221 @@ +#!/usr/bin/env python + +# This script performs the t-test validation for non-bit-for-bit results for the +# CICE model. + +# Written by: Matthew Turner +# Date: October, 2017 + +import netCDF4 as nc +import os +import sys +import numpy as np +import numpy.ma as ma +import itertools +import logging + +def maenumerate(marr): + mask = ~marr.mask.ravel() + for i,m in itertools.izip(np.ndindex(marr.shape[-2:]),mask): + if m: yield i + +import argparse +parser = argparse.ArgumentParser(description='This script performs the T-test for \ + CICE simulations that should be bit-for-bit, but are not.') +parser.add_argument('base_dir',help='Path to the baseline history (iceh_inst*) files. REQUIRED') +parser.add_argument('test_dir',help='Path to the test history (iceh_inst*) files. REQUIRED') +parser.add_argument('-v','--verbose',dest='verbose',help='Print debug output?', \ + action='store_true') + +parser.set_defaults(verbose=False) + +# If no arguments are provided, print the help message +if len(sys.argv)==1: + parser.print_help() + sys.exit(1) + +args = parser.parse_args() + +# Set up the logger +if args.verbose: + logging.basicConfig(level=logging.DEBUG) +else: + logging.basicConfig(level=logging.INFO) +logger = logging.getLogger(__name__) + +# The path to output files for simulation 'a' (the '-bc' simulation) +if args.base_dir.endswith(('history','history/')): + path_a = args.base_dir +else: + path_a = args.base_dir + '/history/' + +# The path to output files for simulation 'b' (the test simulation) +if args.test_dir.endswith(('history','history/')): + path_b = args.test_dir +else: + path_b = args.test_dir + '/history/' + +# Find the number of output files to be read in +files_a = [i for i in os.listdir(path_a+'/') if i.startswith('iceh_inst.')] +files_b = [i for i in os.listdir(path_b+'/') if i.startswith('iceh_inst.')] + +if not len(files_a) == len(files_b): + logger.error("Number of output files for baseline simulation does not match the number" + \ + " of files for the test simulation. Exiting...\n" + \ + "Baseline directory: {}\n".format(path_a) + \ + " # of files: {}\n".format(len(files_a)) + \ + "Test directory: {}\n".format(path_b) + \ + " # of files: {}".format(len(files_b))) + +num_files = len(files_a) +logger.info("Number of files: {}".format(num_files)) + +# Get the array dimensions from the file +nfid = nc.Dataset("{}/{}".format(path_a,files_a[0]),'r') +ni = nfid.dimensions['ni'].size +nj = nfid.dimensions['nj'].size +nfid.close() + +# Pre-allocate a numpy array with a "time" dimension +data_a = np.zeros((num_files,nj,ni)) +data_b = np.zeros((num_files,nj,ni)) + +# Read in the data +var = 'hi' +cnt = 0 +for fname in sorted(files_a): + nfid = nc.Dataset("{}/{}".format(path_a,fname),'r') + fill_value_a = nfid.variables[var]._FillValue + data_a[cnt,:,:] = nfid.variables[var][:] + cnt += 1 + nfid.close() + +cnt = 0 +for fname in sorted(files_b): + nfid = nc.Dataset("{}/{}".format(path_b,fname),'r') + fill_value_b = nfid.variables[var]._FillValue + data_b[cnt,:,:] = nfid.variables[var][:] + cnt += 1 + nfid.close() + +# Calculate the difference and mask the points where the the difference at +# every timestep is 0. +data_d = data_a - data_b +mask_d = np.all(np.equal(data_d,0.),axis=0) +mask_array_a = np.zeros_like(data_d) +mask_array_b = np.zeros_like(data_d) +mask_array_d = np.zeros_like(data_d) +for x,value in np.ndenumerate(mask_d): + i,j = x + mask_array_d[:,i,j] = value + mask_array_a[:,i,j] = value + mask_array_b[:,i,j] = value +data_d = ma.masked_array(data_d,mask=mask_array_d) +data_a = ma.masked_array(data_a,mask=mask_array_a) +data_b = ma.masked_array(data_b,mask=mask_array_b) + +# Calculate the mean of the difference +mean_d = np.mean(data_d,axis=0) +variance_d = np.sum(np.power(data_d - mean_d,2)) / (num_files - 1) + +# Calculate the mean from 1:end-1 and 2:end +mean_nm1_d = np.mean(data_d[:-1,:,:],axis=0) +mean_2n_d = np.mean(data_d[1:,:,:],axis=0) + +# Calculate equation (5) for both simulations +r1_num = np.zeros((nj,ni)) +r1_den1 = np.zeros((nj,ni)) +r1_den2 = np.zeros((nj,ni)) +for i in np.arange(np.size(data_a,axis=0)-1): + r1_num = r1_num + (data_d[i,:,:]-mean_nm1_d[:,:])*(data_d[i+1,:,:]-mean_2n_d[:,:]) + r1_den1 = r1_den1 + np.power(data_d[i,:,:]-mean_nm1_d[:,:],2) + +for i in np.arange(1,np.size(data_a,axis=0)): + r1_den2 = r1_den2 + np.power(data_d[i,:,:] - mean_2n_d[:,:],2) + +r1 = r1_num / np.sqrt(r1_den1*r1_den2) + +# Calculate the effective sample size +n_eff = num_files*( (1.-r1) / (1.+r1) ) +n_eff[n_eff < 2] = 2 +n_eff[n_eff > num_files] = num_files + +# Calculate the t-statistic with n_eff +t_val = mean_d / np.sqrt(variance_d / n_eff) + +# Effective degrees of freedom +df = n_eff - 1 + +# Read in t_crit table +nfid = nc.Dataset("CICE_t_critical_p0.8.nc",'r') +df_table = nfid.variables['df'][:] +t_crit_table = nfid.variables['tcrit'][:] +nfid.close() +t_crit = np.zeros_like(data_d) +for x in maenumerate(data_d): + min_val = np.min(np.abs(df[x]-df_table)) + idx = np.where(np.abs(df[x]-df_table)==min_val) + t_crit[x] = t_crit_table[idx] + +if np.any(abs(t_val) > t_crit): + logger.info('Two-Stage Test Failed') + passed = False + +elif np.all(abs(t_val)<=t_crit) and np.all(n_eff >= 30): + logger.info('Two-Stage Test Passed') + passed = True + +elif np.all(abs(t_val) <= t_crit) and np.any(n_eff < 30): + # Calculate the T-statistic using actual sample size + t_val = mean_d / np.sqrt(variance_d / num_files) + + # Find t_crit from the nearest value on the Lookup Table Test + nfid = nc.Dataset("CICE_t_lookup_p0.8_n1825.nc",'r') + r1_table = nfid.variables['r1'][:] + t_crit_table = nfid.variables['tcrit'][:] + nfid.close() + + for x in maenumerate(data_d): + min_val = np.min(np.abs(r1[x]-r1_table)) + idx = np.where(np.abs(r1[x]-r1_table)==min_val) + t_crit[x] = t_crit_table[idx] + + logger.debug('Max t_crit 1 = {}'.format(np.max(t_crit))) + logger.debug('Min t_crit 1 = {}'.format(np.min(t_crit))) + + if np.any(abs(t_val) > t_crit): + logger.info('Two-Stage Test Failed') + passed = False + elif np.all(abs(t_val) <= t_crit): + logger.info('Two-Stage Test Passed') + passed = True + else: + logger.error('TEST NOT CONCLUSIVE') + +else: + logger.error('TEST NOT CONCLUSIVE') + +# Calculate Taylor Skill Score +#dcorr = np.zeros_like(data_d) +dcorr = np.zeros((nj,ni)) +dcorr[:,:] = 1 +for i,j in maenumerate(data_d): + dcorr[i,j] = np.corrcoef(data_a[:,i,j],data_b[:,i,j])[0,1] +s = np.power( \ + (1+dcorr)*np.std(data_a,axis=0)*np.std(data_b,axis=0) / \ + (np.power(np.std(data_a,axis=0),2) + np.power(np.std(data_b,axis=0),2)) \ + , 2) +s_crit = 0.80 +if np.any(s<0) or np.any(s>1): + logger.error('Skill score out of range') +elif np.all(s > s_crit): + logger.info('Quadratic Skill Test Passed') +else: + logger.info('Quadratic Skill Test Failed') + passed = False + +if passed == False: + sys.exit(1) # exit with an error return code +else: + sys.exit(0) # exit with successfull return code + From cfea740867731a8e2e02f5a5a40752d7eea2b49d Mon Sep 17 00:00:00 2001 From: turner Date: Mon, 20 Nov 2017 21:11:18 +0000 Subject: [PATCH 2/8] Add new option file for the t-test --- configuration/scripts/options/set_nml.ttest | 10 ++++++++++ 1 file changed, 10 insertions(+) create mode 100644 configuration/scripts/options/set_nml.ttest diff --git a/configuration/scripts/options/set_nml.ttest b/configuration/scripts/options/set_nml.ttest new file mode 100644 index 000000000..ba2ae08e7 --- /dev/null +++ b/configuration/scripts/options/set_nml.ttest @@ -0,0 +1,10 @@ +npt = 43800 +dumpfreq = 'm' +dumpfreq_n = 12 +diagfreq = 24 +histfreq = 'd','x','x','x','x' +f_aice = 'd' +f_uvel = 'd' +f_vvel = 'd' +f_hi = 'd' +hist_avg = .false. From 4af2d971ce5b8c8cf24c0096f34cd5803422637f Mon Sep 17 00:00:00 2001 From: turner Date: Wed, 22 Nov 2017 19:15:10 +0000 Subject: [PATCH 3/8] Update t-test script per modifications to the algorithm by Andrew Roberts. Also add the critical t-value lookup table files --- CICE_t_critical_p0.8.nc | Bin 0 -> 15220 bytes CICE_t_lookup_p0.8_n1825.nc | Bin 0 -> 4096 bytes t-test.py | 438 +++++++++++++++++++++--------------- 3 files changed, 254 insertions(+), 184 deletions(-) create mode 100644 CICE_t_critical_p0.8.nc create mode 100644 CICE_t_lookup_p0.8_n1825.nc diff --git a/CICE_t_critical_p0.8.nc b/CICE_t_critical_p0.8.nc new file mode 100644 index 0000000000000000000000000000000000000000..f59775e0a386fcfcc0c0fd05852978fe3629f3d1 GIT binary patch literal 15220 zcma*rb+ld8-7esTyX*0F;_jZnT;!|-5-cPnSO_a}h!B!+a)Qf(28vsuIBa0;2<|Q$ z2=4B5qiu01l$I9mllFe&-ZAb!w_`kg^B!4i&flE-jDBBBI(E(IzyH}#tihVI4qp9* z^_uQSc|7j zpTDGav7cLS=E6e{ZC%j*Kfl*Fc41q*WN~XeCI~mmIqK z$D(zoFKM5%uSq=&T+%l4zYALZGaW|soxHfU zZE@e}ZU6nLzDF!`3uYg*VEUo|*BNb@w$5&A zZC$+jWn=wg*24eJcbzDzU%L92J#O{m*J)pUhX3jRU(fB<@_&xpq8)U6uv4A2aLEi) zSv+@^nJn6&U!P&CpZ$N{)%b15qbTvk#D68L|6QU&HHvnM*qoStBupS<5jm%FE*Dc` z1-Eby6`tZH-lEE1_?DleXmrHdtdFUU-kh!3o_bQgrSZ;l4Hp^iKU#vshrN4oXxqM&jtLEi@B7`xsnpsP_T>@T*nRE#4X&$ z9o)q|+{XhvM1@CL#S=WmGd#x&yu>TK#v8oFJG{pSRQZTc_>8~sIbZTM-|{_w<0pRR zpZtqoqo^gqO>J3|wOAXsqh(#zV|_MYL)`Y3joE}v*^JHEf-TvKt=Wcc*^ce0p_V%8 zX`qoNVw&kdM>^4&E_9_E-RVJ3deNIc^ravD8NfgWF_<9?Wf;R5!AM5213MB(*on~? zMoY>V#TC$f9q7%eUn_@*VlEd{@3F-;?jl_vHuj1NkTUCs~zM`Jwz!ek4DVAIp#BC-M{d zsr*!aCO?yZmVcIik$;hYm4B6=%g^N(@(cN;{8D};zmi|eujM!L8~LsLR(>bHli$nl z~rf0uukKg*xxKjc5;KjlB=FY*`pFZnO|Z~1TetNc~|CV$JL zG>S-A{qIVYuA#4?uc5D@uc@!8uc@!8ucfc0ucfc0udT1GudT1GucNP{ucNP{udA=C zudA=Cucxo4ucxo4udlDKudlDKZ=i3WZ=i3WZ>VpmZ>VpmZ=`ReZ=`ReZ>(>uZ>(>u zZ=!FaZ=!FaZ>n#qZ>n#qZ>DdiZ>DdiZ?12yZ?12yZ=r9YZ=r9YZ>evoZ>evoZ>4Xg zZ>4XgZ>?{wZ>?{wZ=-LcZ=-LcZ>w*sZ>w*sZ>MjkZ>MjkZ?A8!Z?A8!*XT8Rjb5YI z>a}{UUaQyXb$Xp%r`PNCdc9t+H|Py|gWjMw>WzA%-l#X}O?s2wq{n)!$9k+c>&<$z z-mG`fJLnzs4thtuqux>PsCUvk>7DdWdS|_}-dXRgchS4(UGy$`SG}v=Rqv{I)4S>2 z^lo~0y}RCB@2>aId+0s%9(qr`r`}WVsrS-*>Amz`dT+hA-dpdj_tE?4ee^zhU%jv1 zSMRI$)BEZD^nQAOy}#aH@2?Ng2j~Oz0s26FpgvF^s1MQy>4Wq^`e1#qK3E^D57CF{ zL-Zl~P<^OAR3EAj(}(H9^kMpNeYiecAFhwkN9ZH;5&B4dq&`v~sgKe}>7(>f`VRUI z`VRUI`i}aJ`i}aJdeDO&^q?ntq9=Ny@1*af@1*afkJd-)qxI2xi{7HQ=q-Axdq|q< zsXj&@qmR+Y=wtP<`dEFezO%lwzO%lwzKgz#zKgz#K29H}kJHEL3ivW>3iuD z^@;jKeWE@|pQKOHC+U0Zd+U4ad+U?+$@*k{vOYzhqEFGM==ig>Z z>ig;Y>HF#X>HF*Z>-+2b>j&ru=m+Qr=u`En`c!?YexQD!exQD!K24vdPt&LA2k8gt z2k8gtnV#vHp6S!|>H2hix;{gnq0i7~=ri@1`b>SMK1-ja&(dej&!x z>j&$z_1XGteYQSFpQF#w=je0wx%ymvu0Btnr_a;p>4)ft=!fWs==1gY`h0!9eyDz^ zeyDz^zCd4~FVGk03-yKiLVcmWNMEEc(iiE6>4)iu>4)iUdYj&+x9N-Z#rk4>vEHt? z>+O2GzC>T5FVUCihwF#yhwF#yN9afBN9afBN9srFN9srFN9jlDN9jlDN9#xHN9#xH z$LPoC$LPoC$Lh!G$Lh!G$LYuE$LYuE$Lq)I$Lq)IC+H{WC+H{Wxt_-)WaN~nP>s?P zBVrOVa!ORFM(IfrF$oztB`Q>-^mh?42^l#hDpaF%X+%syMox(e)hInVA|@dtr$mKn zl%5h1laP^9qCz!Fe;*N(kdafOLN!WHjfhFe$SF~w8l|U2#3W?ol&DaR($gbi5;Ag1 zRH#Pj84)oF895~?RHO9Fh?sNHOhQIZi3-&yJv$;MAtR?mg=&Hoh)KxEDN&&srRPP&BxK~2s8Eg4^CMysGIB~(s7C1@B4QFUa!ORFM(G6+ zF$oztB`Q>-^umakgp8aL6{=DC$B3AOjGPh`s!@7TL`*_PPKgTDD7`o$CLtrIM1^XU zUJ?KGOhQIZi3-&yEk(p6WaN~nP>s^7BVrOVa!ORFM(H&XF$ozt zB`Q>-^xBA+gp8aL6{=BMM8qUy-bY(s@?B4QFUa!ORFM(NEFF$oztB`Q>-^p=R2gp8aL6{=Bs zYeY;!Mox(e)hN9!A|@dtr$mKnl-?c@laP^9qCz!F?}&&=$jB*Cp&F%kM#Lm!r#u>D>`A2^l#hDpaHNo`{%)jGPh`s!@7xL`*_PPKgTDD7`NtCLtrI zM1^XU-X9T@kdafOLN!Vsh=@tZ$SF~w8l?|L#3W?ol&DaR(uX2q5;Ag1RH#Pj!x1qF z895~?RHL*K5tERSQ=&pON*{@cNyx}4QK1^8k4D5KWaN~nP>s^ZB4QFUa!ORFM(L`E zn1qa+5*4aZ`glZ4LPk!B3e_lmA|fUsBd0`#YLq@15tERSQ=&pON}q~|Nyx}4QK1^8 zPe;TgWaN~nP>s@OB4QFUa!ORFM(MK=F$oztB`Q>-^tp(bgp8aL6{=DCd_+t_Mox(e z)hK--A|@dtr$mKnl)e}dlaP^9qCz!FUy6uH$jB*Cp&F$xN5mv#lViGcP zN>r#u>8lYj2^l#hDpaHNwTPI6jGPh`s!{rSL`*_PPKgTDD19R$CLtrIM1^XUz8Mjd zkdafOLN!X?iikr#u z>4y<92^l#hDpaHNqllPz6md1eOXwQ z)P@KC`gpi!+!x`_`Vrx_XWk08uD5Nt`P#d~O-nxtH(q>dxZ#y~;fC%%hU+i>d$_*# zy>Q(Py~B0mld$sp-@?jbRaiM~SXkLOF0A-utFYpxx5J9#ehVx1io%MnhlCZ;hhh2q zn}y}~-xQW#o`>ZpObE+oejS#N?;n=;{~;`|+d3>?=g6?^yS2lz_pc7io*fyM-TP2j zwz4KHyJW|(?DWR4?3f$FvP1tAmSxw4WqYp`mhJp;ST-sN%ldB}mUY=DEUWjmEw>EI z);}SXe;X6ZKTQkeubvF$Pmc-Z_m2qW*Y6JH7lwuM6F-IW!{3DRy>Tetc1kGU&=kt$ zzlQQvABOU!>xS|LlSBFJ7ee{8LqmD#fKWc6E|ib%8p=!d4&_CqP@cbWD9^qylxL0% zR1dDju4yz`l%-10*xCksP))Ni3Y{F+c6vTrC4^v=Hj z2<2X%gmU)}Lb;3gc3LZxn|p?GqGIc--qI#^F#6T0ipP5QYik$_ftdh?SfEz?R#IY48`ZK zgyOGjh2qa6L-Fa+q4@ZbP<*s$D5^6;@xgjaI@kH-XtolnR9$OKLM-K_b zBLhS6@Nc1b=*3Vxcx@;iI3^VLPYK0+gF|ue_My1j_wM>@DDHeC6n8usireoG#cj8T z;?^5Oaf?0Od_yR1x+N4h-W7@)DxtXk`A}TR}TtBX?iHGIxiGgt_sB!Khrl9m(K~sW!7`4{axA=ic7LkT&!PY?HBb4 z#UJhE!bdTi3(W8j)^Yx?>=lahmWAS6b2xW$D9*W=jY4s@=bdH!XLSw5nWu%~j4$!C zr{9g4o@V~1e#e2K`28~s55*}fs0qc%XM|#@IV^1t#qU1HdQWnuCpo_pUtnA)@+U|` zae_UZAdffi<97+ganCX#6vw`fbsqCE^Fwj;PvltJQ5VxW6i1rtk)uO#ggbD=%upQu zBd219OU!u5m{7F8i#fJ$5Q@duFftTvp53+x&pK>5-oNN0JbzIfiiH(2)(FJ{Gh8q> z6o)#~L(gJBDCWP(adZmBA=Z1yV!U&nz06yHdocGg+@-m`KF1vAETKavX1~k{IJ1MD z>A~lagrfC(3f!euYn}BVhtVYzGwosKxmd>x=QQIE4xwWxru*J>cVznRtQU&R{m+hN zXebWy%!6)X9$iB*?Gwykn!QbPh6ldScn|qt>-WM^8exaD;Ymmq z;~38!YA zF#1b;{zu=+AF%e(`!a$CJY%Pid4iRk#Ul2`ecQ>*lHc$%$&=VeaxU%IXJRc$2h1pV zS9qDbaj(LO%w{*tIQWe2_$#0A4DQN~?!}Jwwc{+t(;q*#gLm#=t~)rZ9qz!H>~IqE z*o)ybWB#MQ#b;oYIgYZ2QRi_KGw{q&1F?@$<~{N&+?kON?;p84_B=9X z!%&QHS4MomQ{0L@jW`AK8R5>37=ynlBWm#M;otB!tMHi}ekmv8=Z7DFa~W<9!_9rT zbr19WVO4y-hMCzgXE^K}tY?_H4RgsqL%eTG59C^-r!d;+rjoU_zKSCXuNOm zK8$4$9odpKLNUlX2btX<=Qrp9ymQcHoQ|^@Gza%^keLh`N*A_eolp$?kxzLIvmAH_ z1sAZCc4o0RDTCeGll1b04q;Ylg!A zjx73HL;pAMjQ*a{|2piU|5+T5bL~GJ=hS};L+FP2_1`!Y{eI;eKEw?BJDPCs|4 z-`V7}(TZ928^8;ZU^^95$w_az?VE>_}MeLbu1iI`1acc<_E*jwKO^Xc0a z>*~7+Q7HQSjX(1a{yy}1ggY^(KAzXdo$B*Dj$j_sn259PGo0R7Umr8?vq30&|C{eH zv)(=fy`RT3dwXVYYwUd`=W_~2vk>>Lw>|gX4eRVZnC>*<9`#-?6uo}IUFzkzy*#&< zwf3_2Ue2S}ja?7JO@GwL=U`|EZ9?p3$(?1a13tsmWJ!XCSAhR;CP zU-^-*Fwd^_+Nt|#H%b!}rFGcen(-rd!^yN+fA1L%Qs z>{^RGcijkQ+2t30;0r#%`@4957w_-#IL@@oo!rDSuEgHE*n5{#G4C!%aTwOn#TvTU zf0s$@Mv61;G8l93(iJo4vOSx#0c(b$vvqW~j?Q0V4xP=v^BcT~``g)>cfOC?FpJIw zSKxEf*;+bVOXm|gn#Ihgm1*q51l)nnX3}{i?p){IbRnjOEisqQYlotfbLsSVd`>(4 zm5(r+POtL<*3`+GI$2YvJGc?^>Ev#7x`aP)CZ})$?tiB?4#B!QSy!ihnaHk;VMm71 zpPqE45wq%KZJn&G(>kH(_!~d-J?7Q%GycTeyvp-DNd@ccXnh@>U&j?(!{uCvxph2^ z-{B5-JOby~(Hc8iW5?oq)8jnCgWQA9MAOYUt0wzvvcIM) zxs(fWUQIr~O{cPyoMW)ZCVOnM$EJB4%ygVvlg~=iWG1j1yI`M93HI4El!5fc`87Ge zrVjY*G}U6SP4?PkuT2}U4r_#>@!vSd#vl0(`)#z}#!s-{M*D4ii`RIG=Xi=$RCoY; zZnWn{XWHmY8&^ZxHHw!jP;H)LJb3`K*x+~CX`{>~44!tGG~CIp+{AS(<65r5JQ~cS!8{ty<19|+6#VUJI01K| z;YivslLj+sn8(4)gQ&a{v)is1}q0A|zRZZw!pgV{7R(}4TYupL{m zIUD1SG?-6=`PBc$zxbIS`HrvofYX)IZNNJi%j_Q~d+n%U#^g z&D_9BmT?VcRew2`a3SY&4rg*2Cvy_!Revl;v4l1jG9UM--aV?Hftl4$Wk1}d`bq4` zZj56r=2mZR^*bV3BA>)4)c*pkiIi1k?q^Q$wzx_|Kxe&Pqd zIr*R5PIg#Tzh9g-*8;fuc>*g_=Sxm=V>&&&zTT>BXxQ^l-n-{mdd;8kAY1)jycYoFjT9>Kh8@8@3b;tp=*Cd|Cn%xjlZ za5YzPIhS$~7hvwS=3aXir*kSN<1W|c9EaK0ntkozEaos4FrT@Yf35k~y4SVlUwa_? zvoDiz2DQ$hb`P9E?KsAgvJ-((IEPy2Q0p9OokMLu-1FL=bf*i>qP7`lQR^&fYuTP{ z*ow{Bl#SSc^;jGCzUDXnjWenF8D~=S1O9&0e8ZQx|258~=2M(YjdQ7SE;aA*HvX2> zyv{2)o0=DJHZ@Q4B&&Fo3J>CZYMf7v^Qmz@HMemyH*!5ISw_J%T*Vb!#wGj_=TzgI zYMfJzbE-Lm(>R5toP@KgISyx4a}-C=P8*AGUNz3EW-iXFrj?mwOk*m}tj3wuII|jO zR% zYEDVXXCSF6uD>VvM702hO)c*n(Q6a9>JXcXa$uw6gjZyz6mx!7zPZ9t% zN>(LLrW6dzDwBDce4-m)Q>_>@}2yh$XI*``3*Ud7l|JLl_8a zlOvu=bM5}Hsu zLf<5rQ76Qz(e79dIkmBj7`I6m663PrE3?Xfi&Uj%-^( zCuA?DS12>-^_#EKn_t|bcSSbQ2l5`$dDc(qV)iv%S<^?~@HC_wwhX4*x`)%xqDIs2 zPB;<04cmxO;Z9=S@e{GO+ee0n9wZ|(j*u}q$H;^``9x%Lib(v5iMOhhOxa#e@M0AS z?zv2+yId!;=HDc9Hryf6g^eV(<33T^wUH&E9fZkvLb#mgB;(F2@`XtcSv~Rt`Kmyd zd>yV&HnbU#jVq1FKWt6NcSp_0_rVrqYvW+DT|0#Av=~ZuXWNpW>EUE=tsU8)Vo!1m z9LV8aqe-sU7;@~2BRQcQPYT{UkyAgokRtIUQc@u%rHkB1*((XD_`4^m8t+A_OMJ+c zIez5YlPRQTEg`iINa_j#Nqu+_x!V>(?ya0gT5Q8e>(LqHQE&w5Xq-hJYa>aQ#cc9C zdoFoNe@0%{Mw6bDg{0RYhV8g>$l<9-E`+*4rs zt`N*aPJ_kfBCu>K2CIoBFeL2^SeKmz8_QC#jXnp%_Lsr%7w5r_l*1?ARDk`hN^lr^ z0Y<5-VD#yWFvhSN#ztNO$K97<-0xRlg7;N$T73_^2lkvE(VtJlO@G>OF(Vi02Ts;|0uq z*bQ?%Uc$Upui&%l*APAI4J?Rz3k!35V9~pG5F7Fbe7?CC{RmybdG` z)`i3cdZ0R>4@x>}P(HNF8CcqY%0%v3f+9-2K%d&vR@kcDI9E9K+e6laQIa*9I`cq zBL#2Z*e||tbmVZzsmOw3hvjhG*B*XZPy_{+hQP09qT%>-V<>DWg2GAHAkTRx9G?9n z6xBb6lIJ(!%&JXr{yTptE@*<&Th_tptr958D}#Kg0hG*ah0-_aP##kY=gmT(G(HH* zVwXeJ&MK&i8x5tOhQg&1gzFxQ;F97TRLNSPdKU|INttl>x*Bd*IKb`E>!F5hh1+^j zP?s7G_2HS&xGxv#dupIQVjr~n{df)FyUHNPH!#NUqLu^ssI}7zs^rF#_95nE6L<8e?G&*X7CVS1%r>*y3a30)7a zN71#NIPoSyF%yTPgU0B+q#eC|E78YxA`-`PoILL>24ps)f5%M3@&zc}L?f0x#(=)1 zD7#>X!L#}>=$SuGQ>5ce(IcF$o`(?)H*kS%7A|bh$LQ=sxNze&jG2EK7gVmsxb0t} zTHTFF-VG?v`4$=d!^j<|!!N9+;R=HTxN=W3emU+Ot~qoBziL{BYnLv@jq?d+zU{$H zhQ|2g$yVHUU4+~BjKHj^XK?rG7~B)pi2DHns&T1JB}lb_-T3 z5HD7az}w%qU~_2)cKPIB?{6EVdYy$*lLw2XR^wFCkrD1vvC2{kW1dMDu56R8uec{Y e;?f$Z)A%>x2Q| num_files] = num_files - -# Calculate the t-statistic with n_eff -t_val = mean_d / np.sqrt(variance_d / n_eff) - -# Effective degrees of freedom -df = n_eff - 1 - -# Read in t_crit table -nfid = nc.Dataset("CICE_t_critical_p0.8.nc",'r') -df_table = nfid.variables['df'][:] -t_crit_table = nfid.variables['tcrit'][:] -nfid.close() -t_crit = np.zeros_like(data_d) -for x in maenumerate(data_d): - min_val = np.min(np.abs(df[x]-df_table)) - idx = np.where(np.abs(df[x]-df_table)==min_val) - t_crit[x] = t_crit_table[idx] - -if np.any(abs(t_val) > t_crit): - logger.info('Two-Stage Test Failed') - passed = False - -elif np.all(abs(t_val)<=t_crit) and np.all(n_eff >= 30): - logger.info('Two-Stage Test Passed') - passed = True - -elif np.all(abs(t_val) <= t_crit) and np.any(n_eff < 30): - # Calculate the T-statistic using actual sample size - t_val = mean_d / np.sqrt(variance_d / num_files) - - # Find t_crit from the nearest value on the Lookup Table Test - nfid = nc.Dataset("CICE_t_lookup_p0.8_n1825.nc",'r') - r1_table = nfid.variables['r1'][:] + + # Pre-allocate a numpy array with a "time" dimension + data_a = np.zeros((num_files,nj,ni)) + data_b = np.zeros((num_files,nj,ni)) + + # Read in the data + var = 'hi' + cnt = 0 + for fname in sorted(files_a): + nfid = nc.Dataset("{}/{}".format(path_a,fname),'r') + fill_value_a = nfid.variables[var]._FillValue + data_a[cnt,:,:] = nfid.variables[var][:] + cnt += 1 + nfid.close() + + cnt = 0 + for fname in sorted(files_b): + nfid = nc.Dataset("{}/{}".format(path_b,fname),'r') + fill_value_b = nfid.variables[var]._FillValue + data_b[cnt,:,:] = nfid.variables[var][:] + cnt += 1 + nfid.close() + + # Calculate the difference and mask the points where the the difference at + # every timestep is 0. + data_d = data_a - data_b + mask_d = np.all(np.equal(data_d,0.),axis=0) + mask_array_a = np.zeros_like(data_d) + mask_array_b = np.zeros_like(data_d) + mask_array_d = np.zeros_like(data_d) + for x,value in np.ndenumerate(mask_d): + i,j = x + mask_array_d[:,i,j] = value + mask_array_a[:,i,j] = value + mask_array_b[:,i,j] = value + data_d = ma.masked_array(data_d,mask=mask_array_d) + data_a = ma.masked_array(data_a,mask=mask_array_a) + data_b = ma.masked_array(data_b,mask=mask_array_b) + + return data_a, data_b, data_d, num_files, path_a, fname + +def two_stage_test(data_a,data_b,num_files,data_d): + # Calculate the mean of the difference + mean_d = np.mean(data_d,axis=0) + variance_d = np.sum(np.power(data_d - mean_d,2)) / (num_files - 1) + + # Calculate the mean from 1:end-1 and 2:end + mean_nm1_d = np.mean(data_d[:-1,:,:],axis=0) + mean_2n_d = np.mean(data_d[1:,:,:],axis=0) + + # Calculate equation (5) for both simulations + r1_num = np.zeros_like(mean_d) + r1_den1 = np.zeros_like(mean_d) + r1_den2 = np.zeros_like(mean_d) + for i in np.arange(np.size(data_a,axis=0)-1): + r1_num = r1_num + (data_d[i,:,:]-mean_nm1_d[:,:])*(data_d[i+1,:,:]-mean_2n_d[:,:]) + r1_den1 = r1_den1 + np.power(data_d[i,:,:]-mean_nm1_d[:,:],2) + + for i in np.arange(1,np.size(data_a,axis=0)): + r1_den2 = r1_den2 + np.power(data_d[i,:,:] - mean_2n_d[:,:],2) + + r1 = r1_num / np.sqrt(r1_den1*r1_den2) + + # Calculate the effective sample size + n_eff = num_files*( (1.-r1) / (1.+r1) ) + n_eff[n_eff < 2] = 2 + n_eff[n_eff > num_files] = num_files + + # Calculate the t-statistic with n_eff + t_val = mean_d / np.sqrt(variance_d / n_eff) + + # Effective degrees of freedom + df = n_eff - 1 + + # Read in t_crit table + nfid = nc.Dataset("CICE_t_critical_p0.8.nc",'r') + df_table = nfid.variables['df'][:] t_crit_table = nfid.variables['tcrit'][:] nfid.close() - + t_crit = np.zeros_like(t_val) for x in maenumerate(data_d): - min_val = np.min(np.abs(r1[x]-r1_table)) - idx = np.where(np.abs(r1[x]-r1_table)==min_val) + min_val = np.min(np.abs(df[x]-df_table)) + idx = np.where(np.abs(df[x]-df_table)==min_val) t_crit[x] = t_crit_table[idx] - - logger.debug('Max t_crit 1 = {}'.format(np.max(t_crit))) - logger.debug('Min t_crit 1 = {}'.format(np.min(t_crit))) if np.any(abs(t_val) > t_crit): logger.info('Two-Stage Test Failed') passed = False - elif np.all(abs(t_val) <= t_crit): + + elif np.all(abs(t_val)<=t_crit) and np.all(n_eff >= 30): logger.info('Two-Stage Test Passed') passed = True + + elif np.all(abs(t_val) <= t_crit) and np.any(n_eff < 30): + # Calculate the T-statistic using actual sample size + t_val = mean_d / np.sqrt(variance_d / num_files) + + # Find t_crit from the nearest value on the Lookup Table Test + nfid = nc.Dataset("CICE_t_lookup_p0.8_n1825.nc",'r') + r1_table = nfid.variables['r1'][:] + t_crit_table = nfid.variables['tcrit'][:] + nfid.close() + + for x in maenumerate(data_d): + min_val = np.min(np.abs(r1[x]-r1_table)) + idx = np.where(np.abs(r1[x]-r1_table)==min_val) + t_crit[x] = t_crit_table[idx] + + if np.any(abs(t_val) > t_crit): + logger.info('Two-Stage Test Failed') + passed = False + elif np.all(abs(t_val) <= t_crit): + logger.info('Two-Stage Test Passed') + passed = True + else: + logger.error('TEST NOT CONCLUSIVE') + passed = False + else: logger.error('TEST NOT CONCLUSIVE') - -else: - logger.error('TEST NOT CONCLUSIVE') + passed = False + return passed # Calculate Taylor Skill Score -#dcorr = np.zeros_like(data_d) -dcorr = np.zeros((nj,ni)) -dcorr[:,:] = 1 -for i,j in maenumerate(data_d): - dcorr[i,j] = np.corrcoef(data_a[:,i,j],data_b[:,i,j])[0,1] -s = np.power( \ - (1+dcorr)*np.std(data_a,axis=0)*np.std(data_b,axis=0) / \ - (np.power(np.std(data_a,axis=0),2) + np.power(np.std(data_b,axis=0),2)) \ - , 2) -s_crit = 0.80 -if np.any(s<0) or np.any(s>1): - logger.error('Skill score out of range') -elif np.all(s > s_crit): - logger.info('Quadratic Skill Test Passed') -else: - logger.info('Quadratic Skill Test Failed') - passed = False - -if passed == False: - sys.exit(1) # exit with an error return code -else: - sys.exit(0) # exit with successfull return code +def skill_test(path_a,fname,data_a,data_b,num_files,tlat): + # First calculate the weight attributed to each grid point (based on Area) + nfid = nc.Dataset("{}/{}".format(path_a,fname),'r') + tarea = nfid.variables['tarea'][:] + nfid.close() + tarea = ma.masked_array(tarea,mask=data_a[0,:,:].mask) + area_weight = tarea / np.sum(tarea) + + weighted_mean_a = 0 + weighted_mean_b = 0 + for i in np.arange(num_files): + weighted_mean_a = weighted_mean_a + np.sum(area_weight*data_a[i,:,:]) + weighted_mean_b = weighted_mean_b + np.sum(area_weight*data_b[i,:,:]) + + weighted_mean_a = weighted_mean_a / num_files + weighted_mean_b = weighted_mean_b / num_files + + nonzero_weights = np.count_nonzero(area_weight) + area_var_a = 0 + area_var_b = 0 + for t in np.arange(num_files): + area_var_a = area_var_a + np.sum(area_weight*np.power(data_a[t,:,:]-weighted_mean_a,2)) + area_var_b = area_var_b + np.sum(area_weight*np.power(data_b[t,:,:]-weighted_mean_b,2)) + + area_var_a = nonzero_weights / (num_files * nonzero_weights - 1.) * area_var_a + area_var_b = nonzero_weights / (num_files * nonzero_weights - 1.) * area_var_b + std_a = np.sqrt(area_var_a) + std_b = np.sqrt(area_var_b) + + combined_cov = 0 + for i in np.arange(num_files): + combined_cov = combined_cov + np.sum(area_weight*(data_a[i,:,:]-weighted_mean_a)*\ + (data_b[i,:,:]-weighted_mean_b)) + + combined_cov = nonzero_weights / (num_files * nonzero_weights - 1.) * combined_cov + + weighted_r = combined_cov / (std_a*std_b) + + s = np.power((1+weighted_r)*(std_a*std_b)/\ + (area_var_a + area_var_b),2) + + logger.debug('s = {}'.format(s)) + + s_crit = 0.99 + if s<0 or s>1: + logger.error('Skill score out of range') + passed = False + elif s > s_crit: + logger.info('Quadratic Skill Test Passed') + passed = True + else: + logger.info('Quadratic Skill Test Failed') + passed = False + return passed + +if __name__ == "__main__": + import argparse + parser = argparse.ArgumentParser(description='This script performs the T-test for \ + CICE simulations that should be bit-for-bit, but are not.') + parser.add_argument('base_dir',help='Path to the baseline history (iceh_inst*) files. REQUIRED') + parser.add_argument('test_dir',help='Path to the test history (iceh_inst*) files. REQUIRED') + parser.add_argument('-v','--verbose',dest='verbose',help='Print debug output?', \ + action='store_true') + + parser.set_defaults(verbose=False) + + # If no arguments are provided, print the help message + if len(sys.argv)==1: + parser.print_help() + sys.exit(1) + + args = parser.parse_args() + + # Set up the logger + if args.verbose: + logging.basicConfig(level=logging.DEBUG) + else: + logging.basicConfig(level=logging.INFO) + logger = logging.getLogger(__name__) + + data_a, data_b, data_d, num_files, path_a, fname = read_data(args.base_dir,args.test_dir) + # Run the two-stage test + passed = two_stage_test(data_a,data_b,num_files,data_d) + + nfid = nc.Dataset("{}/{}".format(path_a,fname),'r') + tlat = nfid.variables['TLAT'][:] + nfid.close() + + mask_tlat = tlat < 0 + mask_nh = np.zeros_like(data_a) + mask_sh = np.zeros_like(data_a) + for (i,j),value in np.ndenumerate(mask_tlat): + mask_nh[:,i,j] = value + mask_sh[:,i,j] = not value + + # Run skill test on northern hemisphere + data_nh_a = ma.masked_array(data_a,mask=mask_nh) + data_nh_b = ma.masked_array(data_b,mask=mask_nh) + passed_nh = skill_test(path_a,fname,data_nh_a,data_nh_b,num_files,tlat) + + # Run skill test on southern hemisphere + data_sh_a = ma.masked_array(data_a,mask=mask_sh) + data_sh_b = ma.masked_array(data_b,mask=mask_sh) + passed_sh = skill_test(path_a,fname,data_sh_a,data_sh_b,num_files,tlat) + + passed_skill = passed_nh and passed_sh + + if not passed and not passed_skill: + sys.exit(1) # exit with an error return code + else: + sys.exit(0) # exit with successfull return code + From b945eecf574b5f2eff5672a0cacb4917944bbb9a Mon Sep 17 00:00:00 2001 From: turner Date: Wed, 22 Nov 2017 13:30:24 -0600 Subject: [PATCH 4/8] move QC files to a new QC folder within configuration/scripts/tests --- .../scripts/tests/QC/CICE_t_critical_p0.8.nc | Bin .../scripts/tests/QC/CICE_t_lookup_p0.8_n1825.nc | Bin .../scripts/tests/QC/cice.t-test.py | 15 +++++++++------ 3 files changed, 9 insertions(+), 6 deletions(-) rename CICE_t_critical_p0.8.nc => configuration/scripts/tests/QC/CICE_t_critical_p0.8.nc (100%) rename CICE_t_lookup_p0.8_n1825.nc => configuration/scripts/tests/QC/CICE_t_lookup_p0.8_n1825.nc (100%) rename t-test.py => configuration/scripts/tests/QC/cice.t-test.py (95%) diff --git a/CICE_t_critical_p0.8.nc b/configuration/scripts/tests/QC/CICE_t_critical_p0.8.nc similarity index 100% rename from CICE_t_critical_p0.8.nc rename to configuration/scripts/tests/QC/CICE_t_critical_p0.8.nc diff --git a/CICE_t_lookup_p0.8_n1825.nc b/configuration/scripts/tests/QC/CICE_t_lookup_p0.8_n1825.nc similarity index 100% rename from CICE_t_lookup_p0.8_n1825.nc rename to configuration/scripts/tests/QC/CICE_t_lookup_p0.8_n1825.nc diff --git a/t-test.py b/configuration/scripts/tests/QC/cice.t-test.py similarity index 95% rename from t-test.py rename to configuration/scripts/tests/QC/cice.t-test.py index b79db1526..c0fac089e 100755 --- a/t-test.py +++ b/configuration/scripts/tests/QC/cice.t-test.py @@ -176,7 +176,7 @@ def two_stage_test(data_a,data_b,num_files,data_d): return passed # Calculate Taylor Skill Score -def skill_test(path_a,fname,data_a,data_b,num_files,tlat): +def skill_test(path_a,fname,data_a,data_b,num_files,tlat,hemisphere): # First calculate the weight attributed to each grid point (based on Area) nfid = nc.Dataset("{}/{}".format(path_a,fname),'r') tarea = nfid.variables['tarea'][:] @@ -221,13 +221,13 @@ def skill_test(path_a,fname,data_a,data_b,num_files,tlat): s_crit = 0.99 if s<0 or s>1: - logger.error('Skill score out of range') + logger.error('Skill score out of range for {} Hemisphere'.format(hemisphere)) passed = False elif s > s_crit: - logger.info('Quadratic Skill Test Passed') + logger.info('Quadratic Skill Test Passed for {} Hemisphere'.format(hemisphere)) passed = True else: - logger.info('Quadratic Skill Test Failed') + logger.info('Quadratic Skill Test Failed for {} Hemisphere'.format(hemisphere)) passed = False return passed @@ -275,17 +275,20 @@ def skill_test(path_a,fname,data_a,data_b,num_files,tlat): # Run skill test on northern hemisphere data_nh_a = ma.masked_array(data_a,mask=mask_nh) data_nh_b = ma.masked_array(data_b,mask=mask_nh) - passed_nh = skill_test(path_a,fname,data_nh_a,data_nh_b,num_files,tlat) + passed_nh = skill_test(path_a,fname,data_nh_a,data_nh_b,num_files,tlat,'Northern') # Run skill test on southern hemisphere data_sh_a = ma.masked_array(data_a,mask=mask_sh) data_sh_b = ma.masked_array(data_b,mask=mask_sh) - passed_sh = skill_test(path_a,fname,data_sh_a,data_sh_b,num_files,tlat) + passed_sh = skill_test(path_a,fname,data_sh_a,data_sh_b,num_files,tlat,'Southern') passed_skill = passed_nh and passed_sh + logger.info('') if not passed and not passed_skill: + logger.error('Quality Control Test FAILED') sys.exit(1) # exit with an error return code else: + logger.info('Quality Control Test PASSED') sys.exit(0) # exit with successfull return code From d721fc788c7ab41d3384da0fe8f1a91e00a0083b Mon Sep 17 00:00:00 2001 From: turner Date: Wed, 22 Nov 2017 13:32:14 -0600 Subject: [PATCH 5/8] Update t-test script per the new location of the lookup table files --- configuration/scripts/tests/QC/cice.t-test.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/configuration/scripts/tests/QC/cice.t-test.py b/configuration/scripts/tests/QC/cice.t-test.py index c0fac089e..449abf17a 100755 --- a/configuration/scripts/tests/QC/cice.t-test.py +++ b/configuration/scripts/tests/QC/cice.t-test.py @@ -127,7 +127,7 @@ def two_stage_test(data_a,data_b,num_files,data_d): df = n_eff - 1 # Read in t_crit table - nfid = nc.Dataset("CICE_t_critical_p0.8.nc",'r') + nfid = nc.Dataset("configuration/scripts/tests/QC/CICE_t_critical_p0.8.nc",'r') df_table = nfid.variables['df'][:] t_crit_table = nfid.variables['tcrit'][:] nfid.close() @@ -150,7 +150,7 @@ def two_stage_test(data_a,data_b,num_files,data_d): t_val = mean_d / np.sqrt(variance_d / num_files) # Find t_crit from the nearest value on the Lookup Table Test - nfid = nc.Dataset("CICE_t_lookup_p0.8_n1825.nc",'r') + nfid = nc.Dataset("configuration/scripts/tests/QC/CICE_t_lookup_p0.8_n1825.nc",'r') r1_table = nfid.variables['r1'][:] t_crit_table = nfid.variables['tcrit'][:] nfid.close() From 8243a68815a45b3fc11fbc9d835776c857de4b37 Mon Sep 17 00:00:00 2001 From: turner Date: Wed, 22 Nov 2017 14:34:40 -0600 Subject: [PATCH 6/8] Update cice.t-test.py to be compatible with Python 3 in addition to Python 2 --- configuration/scripts/tests/QC/cice.t-test.py | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/configuration/scripts/tests/QC/cice.t-test.py b/configuration/scripts/tests/QC/cice.t-test.py index 449abf17a..c1a0d9e9d 100755 --- a/configuration/scripts/tests/QC/cice.t-test.py +++ b/configuration/scripts/tests/QC/cice.t-test.py @@ -11,13 +11,17 @@ import sys import numpy as np import numpy.ma as ma -import itertools import logging def maenumerate(marr): mask = ~marr.mask.ravel() - for i,m in itertools.izip(np.ndindex(marr.shape[-2:]),mask): - if m: yield i + try: # Python 2 + import itertools + for i,m in itertools.izip(np.ndindex(marr.shape[-2:]),mask): + if m: yield i + except: # Python 3 + for i,m in zip(np.ndindex(marr.shape[-2:]),mask): + if m: yield i def read_data(base_dir,test_dir): # The path to output files for simulation 'a' (the '-bc' simulation) From 7d63f44b95c5c9cdc8be07e195e7f051b9887e95 Mon Sep 17 00:00:00 2001 From: turner Date: Wed, 22 Nov 2017 14:35:04 -0600 Subject: [PATCH 7/8] Update documentation with the practical compliance procedure for the t-test --- doc/source/cice_3_user_guide.rst | 29 +++++++++++++++++++++++++++-- 1 file changed, 27 insertions(+), 2 deletions(-) diff --git a/doc/source/cice_3_user_guide.rst b/doc/source/cice_3_user_guide.rst index 1c074a515..de6356b5b 100644 --- a/doc/source/cice_3_user_guide.rst +++ b/doc/source/cice_3_user_guide.rst @@ -2184,8 +2184,33 @@ the first and last of these values corresponding to :math:`\sigma` and Practical Testing Procedure *************************** -To be placed here: Write up of how to actually do this test within the -testing software to be added by Elizabeth, Rick, Matt, Tony et al.... +The CICE code compliance test is performed by running a python script (cice.t-test.py). +In order to run the script, the following requirements must be met: +- Python v2.7 or later +- netCDF Python package +- numpy Python package + +In order to generate the files necessary for the compliance test, test cases should be +created with the ``ttest`` option (i.e., ``-s ttest``) when running create.case. This +option results in daily, non-averaged history files being written for a 5 year simulation. + +To run the compliance test: +.. code-block:: bash + + cp configuration/scripts/tests/QC/cice.t-test.py . + ./cice.t-test.py /path/to/baseline/history /path/to/test/history + +The script will produce output similar to: + + INFO:__main__:Number of files: 1825 + INFO:__main__:Two-Stage Test Passed + INFO:__main__:Quadratic Skill Test Passed for Northern Hemisphere + INFO:__main__:Quadratic Skill Test Passed for Southern Hemisphere + INFO:__main__: + INFO:__main__:Quality Control Test PASSED + +Additionally, the exit code from the test (``echo $?``) will be 0 if the test passed, +and 1 if the test failed. Implementation notes: 1) Provide a pass/fail on each of the confidence intervals, 2) Facilitate output of a bitmap for each test so that From 1ed832f757b6707d780b1cf4fc437820828e920b Mon Sep 17 00:00:00 2001 From: "Matthew Turner (Contractor)" Date: Wed, 22 Nov 2017 14:48:29 -0600 Subject: [PATCH 8/8] Update formatting for compliance test practical testing procedure documentation --- doc/source/cice_3_user_guide.rst | 20 +++++++++++--------- 1 file changed, 11 insertions(+), 9 deletions(-) diff --git a/doc/source/cice_3_user_guide.rst b/doc/source/cice_3_user_guide.rst index de6356b5b..d28fd3874 100644 --- a/doc/source/cice_3_user_guide.rst +++ b/doc/source/cice_3_user_guide.rst @@ -2186,15 +2186,17 @@ Practical Testing Procedure The CICE code compliance test is performed by running a python script (cice.t-test.py). In order to run the script, the following requirements must be met: -- Python v2.7 or later -- netCDF Python package -- numpy Python package + +* Python v2.7 or later +* netCDF Python package +* numpy Python package In order to generate the files necessary for the compliance test, test cases should be created with the ``ttest`` option (i.e., ``-s ttest``) when running create.case. This option results in daily, non-averaged history files being written for a 5 year simulation. To run the compliance test: + .. code-block:: bash cp configuration/scripts/tests/QC/cice.t-test.py . @@ -2202,12 +2204,12 @@ To run the compliance test: The script will produce output similar to: - INFO:__main__:Number of files: 1825 - INFO:__main__:Two-Stage Test Passed - INFO:__main__:Quadratic Skill Test Passed for Northern Hemisphere - INFO:__main__:Quadratic Skill Test Passed for Southern Hemisphere - INFO:__main__: - INFO:__main__:Quality Control Test PASSED + | \INFO:__main__:Number of files: 1825 + | \INFO:__main__:Two-Stage Test Passed + | \INFO:__main__:Quadratic Skill Test Passed for Northern Hemisphere + | \INFO:__main__:Quadratic Skill Test Passed for Southern Hemisphere + | \INFO:__main__: + | \INFO:__main__:Quality Control Test PASSED Additionally, the exit code from the test (``echo $?``) will be 0 if the test passed, and 1 if the test failed.