Skip to content

Update py-bkerror Python scripts from Python 2 to Python 3 #93

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
Copilot wants to merge 3 commits into
base: develop
Choose a base branch
from
Draft

Update py-bkerror Python scripts from Python 2 to Python 3 #93

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
a814a53
Update py-bkerror scripts from Python 2 to Python 3
Copilot Apr 15, 2026
0846456
Improve _interp2d docstring to document kind parameter
Copilot Apr 15, 2026
a6d7f14
Fix compile.csh: use --fcompiler=$compiler consistently; add ifx support
Copilot Apr 15, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions src/NMC_Bkerror/py-bkerror/avg_bkerror.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,11 +4,11 @@
nmon = 0
for mon in ['jan','feb','mar','apr','may','june','july','aug','sept','oct','nov','dec']:
filename = '../../../../../staticB/24h/global_berror.l127y770.f77_%ssmooth0p5' % (mon,)
print nmon,filename
print(nmon, filename)
nsig,nlat,nlon = bkerror.get_header(filename)
ivar,agvin,bgvin,wgvin,corzin,hscalesin,vscalesin,corq2in,corsstin,hsstin,corpin,hscalespin = bkerror.get_bkerror(filename,nsig,nlat,nlon)
if not nmon:
print 'initalize arrays'
print('initalize arrays')
agvout = np.zeros(agvin.shape, agvin.dtype)
bgvout = np.zeros(bgvin.shape, bgvin.dtype)
wgvout = np.zeros(wgvin.shape, wgvin.dtype)
Expand All @@ -34,7 +34,7 @@
nmon += 1

filename = '../../../../../staticB/24h/global_berror.l127y770.f77_annmeansmooth0p5'
print filename
print(filename)
bkerror.put_bkerror(filename,ivar,\
agvout,bgvout,wgvout,\
corzout,hscalesout,vscalesout,\
Expand Down
18 changes: 14 additions & 4 deletions src/NMC_Bkerror/py-bkerror/compile.csh
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,9 @@ else if ( $compiler == "intel" ) then
else if ( $compiler == "intelem" ) then
set FC = "ifort"
set FFLAGS = "-g -traceback"
else if ( $compiler == "ifx" ) then
set FC = "ifx"
set FFLAGS = "-g -traceback"
else
echo "UNKNOWN compiler = $compiler"
echo "USING GNU95 compiler instead"
Expand All @@ -39,10 +42,17 @@ if ( $todo == "clean" ) then
else if ( $todo == "build" ) then

$FC -c $FFLAGS str2arr2str.f90
f2py -c splat.F -m splat --fcompiler=$compiler
f2py -c str2arr2str.f90 -m str2arr2str --fcompiler=$compiler
f2py -m bkerror -h bkerror.pyf bkerror.f90
f2py -c --fcompiler=$FC bkerror.pyf bkerror.f90 str2arr2str.o
if ( $compiler == "ifx" ) then
f2py -c splat.F -m splat --f90exec=ifx --f77exec=ifx
f2py -c str2arr2str.f90 -m str2arr2str --f90exec=ifx --f77exec=ifx
f2py -m bkerror -h bkerror.pyf bkerror.f90
f2py -c --f90exec=ifx --f77exec=ifx bkerror.pyf bkerror.f90 str2arr2str.o
else
f2py -c splat.F -m splat --fcompiler=$compiler
f2py -c str2arr2str.f90 -m str2arr2str --fcompiler=$compiler
f2py -m bkerror -h bkerror.pyf bkerror.f90
f2py -c --fcompiler=$compiler bkerror.pyf bkerror.f90 str2arr2str.o
endif

else if ( $todo == "test" ) then

Expand Down
94 changes: 47 additions & 47 deletions src/NMC_Bkerror/py-bkerror/error_bkerror.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
#!/usr/bin/env python
#!/usr/bin/env python3

import sys
import numpy as np
Expand All @@ -9,59 +9,59 @@

nsig1,nlat1,nlon1 = bkerror.get_header(f1name)
ivar1,agvin1,bgvin1,wgvin1,corzin1,hscalesin1,vscalesin1,corq2in1,corsstin1,hsstin1,corpin1,hscalespin1 = bkerror.get_bkerror(f1name,nsig1,nlat1,nlon1)
var1 = (ivar1.tostring()).replace('\x00','')[:-1].split('|')
var1 = ivar1.tobytes().decode('ascii').replace('\x00','')[:-1].split('|')

nsig2,nlat2,nlon2 = bkerror.get_header(f2name)
ivar2,agvin2,bgvin2,wgvin2,corzin2,hscalesin2,vscalesin2,corq2in2,corsstin2,hsstin2,corpin2,hscalespin2 = bkerror.get_bkerror(f2name,nsig2,nlat2,nlon2)
var2 = (ivar2.tostring()).replace('\x00','')[:-1].split('|')
var2 = ivar2.tobytes().decode('ascii').replace('\x00','')[:-1].split('|')

# Print some info about the file
print 'info from %s' % f1name
print 'nsig = %d, nlat = %d, nlon = %d, nvar = %d' % (nsig1,nlat1,nlon1,len(var1))
print 'variables in %s' % f1name
print ', '.join(var1)
print 'agvin.shape: ', agvin1.shape
print 'bgvin.shape: ', bgvin1.shape
print 'wgvin.shape: ', wgvin1.shape
print 'corzin.shape: ', corzin1.shape
print 'hscalesin.shape: ', hscalesin1.shape
print 'vscalesin.shape: ', vscalesin1.shape
print 'corq2in.shape: ', corq2in1.shape
print 'corsstin.shape: ', corsstin1.shape
print 'hsstin.shape: ', hsstin1.shape
print 'corpin.shape: ', corpin1.shape
print 'hscalespin.shape: ', hscalespin1.shape
print('info from %s' % f1name)
print('nsig = %d, nlat = %d, nlon = %d, nvar = %d' % (nsig1,nlat1,nlon1,len(var1)))
print('variables in %s' % f1name)
print(', '.join(var1))
print('agvin.shape: ', agvin1.shape)
print('bgvin.shape: ', bgvin1.shape)
print('wgvin.shape: ', wgvin1.shape)
print('corzin.shape: ', corzin1.shape)
print('hscalesin.shape: ', hscalesin1.shape)
print('vscalesin.shape: ', vscalesin1.shape)
print('corq2in.shape: ', corq2in1.shape)
print('corsstin.shape: ', corsstin1.shape)
print('hsstin.shape: ', hsstin1.shape)
print('corpin.shape: ', corpin1.shape)
print('hscalespin.shape: ', hscalespin1.shape)

print
print()

# Print some info about the file
print 'info from %s' % f2name
print 'nsig = %d, nlat = %d, nlon = %d, nvar = %d' % (nsig2,nlat2,nlon2,len(var2))
print 'variables in %s' % f2name
print ', '.join(var2)
print 'agvin.shape: ', agvin2.shape
print 'bgvin.shape: ', bgvin2.shape
print 'wgvin.shape: ', wgvin2.shape
print 'corzin.shape: ', corzin2.shape
print 'hscalesin.shape: ', hscalesin2.shape
print 'vscalesin.shape: ', vscalesin2.shape
print 'corq2in.shape: ', corq2in2.shape
print 'corsstin.shape: ', corsstin2.shape
print 'hsstin.shape: ', hsstin2.shape
print 'corpin.shape: ', corpin2.shape
print 'hscalespin.shape: ', hscalespin2.shape
print('info from %s' % f2name)
print('nsig = %d, nlat = %d, nlon = %d, nvar = %d' % (nsig2,nlat2,nlon2,len(var2)))
print('variables in %s' % f2name)
print(', '.join(var2))
print('agvin.shape: ', agvin2.shape)
print('bgvin.shape: ', bgvin2.shape)
print('wgvin.shape: ', wgvin2.shape)
print('corzin.shape: ', corzin2.shape)
print('hscalesin.shape: ', hscalesin2.shape)
print('vscalesin.shape: ', vscalesin2.shape)
print('corq2in.shape: ', corq2in2.shape)
print('corsstin.shape: ', corsstin2.shape)
print('hsstin.shape: ', hsstin2.shape)
print('corpin.shape: ', corpin2.shape)
print('hscalespin.shape: ', hscalespin2.shape)

print
print 'differences'
print 'agvin ', np.abs(agvin1-agvin2).max()
print 'bgvin ', np.abs(bgvin1-bgvin2).max()
print 'wgvin ', np.abs(wgvin1-wgvin2).max()
print 'corzin ', np.abs(corzin1-corzin2).max()
print 'hscalesin ', np.abs(hscalesin1-hscalesin2).max()
print 'vscalesin ', np.abs(vscalesin1-vscalesin2).max()
print 'corq2in ', np.abs(corq2in1-corq2in2).max()
print 'corsstin ', np.abs(corsstin1-corsstin2).max()
print 'hsstin ', np.abs(hsstin1-hsstin2).max()
print 'corpin ', np.abs(corpin1-corpin2).max()
print 'hscalespin ', np.abs(hscalespin1-hscalespin2).max()
print()
print('differences')
print('agvin ', np.abs(agvin1-agvin2).max())
print('bgvin ', np.abs(bgvin1-bgvin2).max())
print('wgvin ', np.abs(wgvin1-wgvin2).max())
print('corzin ', np.abs(corzin1-corzin2).max())
print('hscalesin ', np.abs(hscalesin1-hscalesin2).max())
print('vscalesin ', np.abs(vscalesin1-vscalesin2).max())
print('corq2in ', np.abs(corq2in1-corq2in2).max())
print('corsstin ', np.abs(corsstin1-corsstin2).max())
print('hsstin ', np.abs(hsstin1-hsstin2).max())
print('corpin ', np.abs(corpin1-corpin2).max())
print('hscalespin ', np.abs(hscalespin1-hscalespin2).max())

87 changes: 51 additions & 36 deletions src/NMC_Bkerror/py-bkerror/interp_bkerror.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,8 +1,8 @@
#!/usr/bin/env python
#!/usr/bin/env python3

import sys
import numpy as np
from scipy.interpolate import interp1d, interp2d
from scipy.interpolate import interp1d, RegularGridInterpolator
from argparse import ArgumentParser, ArgumentDefaultsHelpFormatter
from bkerror import bkerror
from splat import splat
Expand All @@ -18,7 +18,7 @@ def __init__(self,filename):
'''
nsig,nlat,nlon = bkerror.get_header(filename)
ivar,agvin,bgvin,wgvin,corzin,hscalesin,vscalesin,corq2in,corsstin,hsstin,corpin,hscalespin = bkerror.get_bkerror(filename,nsig,nlat,nlon)
var = (ivar.tostring()).replace('\x00','')[:-1].split('|')
var = ivar.tobytes().decode('ascii').replace('\x00','')[:-1].split('|')

self.filename = filename

Expand Down Expand Up @@ -49,22 +49,22 @@ def print_summary(self):
Print a summary of the GSI background error file
'''

print
print 'file = %s' % self.filename
print 'nsig = %d, nlat = %d, nlon = %d, nvar = %d' % (self.nsig,self.nlat,self.nlon,len(self.var))
print 'variables = %s' % ', '.join(self.var)
print 'agv.shape: ', self.agvin.shape
print 'bgv.shape: ', self.bgvin.shape
print 'wgv.shape: ', self.wgvin.shape
print 'corz.shape: ', self.corzin.shape
print 'hscales.shape: ', self.hscalesin.shape
print 'vscales.shape: ', self.vscalesin.shape
print 'corq2.shape: ', self.corq2in.shape
print 'corsst.shape: ', self.corsstin.shape
print 'hsst.shape: ', self.hsstin.shape
print 'corp.shape: ', self.corpin.shape
print 'hscalesp.shape: ', self.hscalespin.shape
print
print()
print('file = %s' % self.filename)
print('nsig = %d, nlat = %d, nlon = %d, nvar = %d' % (self.nsig,self.nlat,self.nlon,len(self.var)))
print('variables = %s' % ', '.join(self.var))
print('agv.shape: ', self.agvin.shape)
print('bgv.shape: ', self.bgvin.shape)
print('wgv.shape: ', self.wgvin.shape)
print('corz.shape: ', self.corzin.shape)
print('hscales.shape: ', self.hscalesin.shape)
print('vscales.shape: ', self.vscalesin.shape)
print('corq2.shape: ', self.corq2in.shape)
print('corsst.shape: ', self.corsstin.shape)
print('hsst.shape: ', self.hsstin.shape)
print('corp.shape: ', self.corpin.shape)
print('hscalesp.shape: ', self.hscalespin.shape)
print()

return

Expand Down Expand Up @@ -107,53 +107,68 @@ def print_summary(self):
slon_n = np.linspace(0.,360.,gsi_n.nlon,endpoint=False)
slat_n = 180. / np.arccos(-1.) * np.arcsin(glat_n[::-1])

print 'Interpolate from %d to %d' % (gsi.nlat, gsi_n.nlat)
print('Interpolate from %d to %d' % (gsi.nlat, gsi_n.nlat))

def _interp2d(x, y, z, kind='linear'):
'''2D interpolation replacing deprecated scipy.interpolate.interp2d.
x: 1D array of column coordinates (must be strictly ascending)
y: 1D array of row coordinates (must be strictly ascending)
z: 2D array of shape (len(y), len(x))
kind: interpolation method - 'linear' or 'cubic' (requires scipy >= 1.9)
Returns a callable f(x_new, y_new) producing shape (len(y_new), len(x_new)).
'''
f = RegularGridInterpolator((y, x), z, method=kind, bounds_error=False, fill_value=None)
def wrapper(x_new, y_new):
yy, xx = np.meshgrid(y_new, x_new, indexing='ij')
pts = np.column_stack([yy.ravel(), xx.ravel()])
return f(pts).reshape(len(y_new), len(x_new))
return wrapper

tmp = gsi.agvin.reshape(gsi.nlat,-1)
f = interp2d(ssig2,slat,tmp,kind=interp_kind)
f = _interp2d(ssig2,slat,tmp,kind=interp_kind)
tmp_n = f(ssig2,slat_n)
gsi_n.agvin = np.array(tmp_n.reshape(gsi_n.nlat,gsi_n.nsig,gsi_n.nsig),dtype=np.float32)

f = interp2d(ssig,slat,gsi.bgvin,kind=interp_kind)
f = _interp2d(ssig,slat,gsi.bgvin,kind=interp_kind)
tmp_n = f(ssig,slat_n)
gsi_n.bgvin = np.array(tmp_n,dtype=np.float32)

f = interp2d(ssig,slat,gsi.wgvin,kind=interp_kind)
f = _interp2d(ssig,slat,gsi.wgvin,kind=interp_kind)
tmp_n = f(ssig,slat_n)
gsi_n.wgvin = np.array(tmp_n,dtype=np.float32)

tmp = gsi.corzin.reshape(gsi.nlat,-1)
f = interp2d(ssig3,slat,tmp,kind=interp_kind)
f = _interp2d(ssig3,slat,tmp,kind=interp_kind)
tmp_n = f(ssig3,slat_n)
gsi_n.corzin = np.array(tmp_n.reshape(gsi_n.nlat,gsi_n.nsig,6),dtype=np.float32)

tmp = gsi.hscalesin.reshape(gsi.nlat,-1)
f = interp2d(ssig3,slat,tmp,kind=interp_kind)
f = _interp2d(ssig3,slat,tmp,kind=interp_kind)
tmp_n = f(ssig3,slat_n)
gsi_n.hscalesin = np.array(tmp_n.reshape(gsi_n.nlat,gsi_n.nsig,6),dtype=np.float32)

tmp = gsi.vscalesin.reshape(gsi.nlat,-1)
f = interp2d(ssig3,slat,tmp,kind=interp_kind)
f = _interp2d(ssig3,slat,tmp,kind=interp_kind)
tmp_n = f(ssig3,slat_n)
gsi_n.vscalesin = np.array(tmp_n.reshape(gsi_n.nlat,gsi_n.nsig,6),dtype=np.float32)

f = interp2d(ssig,slat,gsi.corq2in,kind=interp_kind)
f = _interp2d(ssig,slat,gsi.corq2in,kind=interp_kind)
tmp_n = f(ssig,slat_n)
gsi_n.corq2in = np.array(tmp_n,dtype=np.float32)

f = interp2d(slon,slat,gsi.corsstin,kind=interp_kind)
f = _interp2d(slon,slat,gsi.corsstin,kind=interp_kind)
tmp_n = f(slon_n,slat_n)
gsi_n.corsstin = np.array(tmp_n,dtype=np.float32)

f = interp2d(slon,slat,gsi.hsstin,kind=interp_kind)
f = _interp2d(slon,slat,gsi.hsstin,kind=interp_kind)
tmp_n = f(slon_n,slat_n)
gsi_n.hsstin = np.array(tmp_n,dtype=np.float32)

f = interp1d(slat.astype(np.float),gsi.corpin.astype(np.float),kind=interp_kind,fill_value="extrapolate")
f = interp1d(slat.astype(float),gsi.corpin.astype(float),kind=interp_kind,fill_value="extrapolate")
tmp_n = f(slat_n)
gsi_n.corpin = np.array(tmp_n,dtype=np.float32)

f = interp1d(slat.astype(np.float),gsi.hscalespin.astype(np.float),kind=interp_kind,fill_value="extrapolate")
f = interp1d(slat.astype(float),gsi.hscalespin.astype(float),kind=interp_kind,fill_value="extrapolate")
tmp_n = f(slat_n)
gsi_n.hscalespin = np.array(tmp_n,dtype=np.float32)

Expand All @@ -169,8 +184,8 @@ def print_summary(self):
gsi_rn = GSIbkgerr(gsi_n.filename)
gsi_rn.print_summary()

print 'differences'
print np.abs(gsi_n.agvin-gsi_rn.agvin).max()
print np.abs(gsi_n.bgvin-gsi_rn.bgvin).max()
print np.abs(gsi_n.wgvin-gsi_rn.wgvin).max()
print np.abs(gsi_n.hscalesin-gsi_rn.hscalesin).max()
print('differences')
print(np.abs(gsi_n.agvin-gsi_rn.agvin).max())
print(np.abs(gsi_n.bgvin-gsi_rn.bgvin).max())
print(np.abs(gsi_n.wgvin-gsi_rn.wgvin).max())
print(np.abs(gsi_n.hscalesin-gsi_rn.hscalesin).max())
Loading