Update to README on how to setup PYTHONPATH and what scripts use depr…

…ecated units package. Small fixes to other python scripts so they run.

helpers/README.md

@@ -1,15 +1,28 @@
 # Install Python Dependencies
 The following instructions and dependecy files work for the core ICAR scripts.
-Tools in `make_domain.py` and ccsm, cesm, cmip, erai, and wrf directories will require Bunch and mygis packages as well.
+Tools in `make_domain.py` and ccsm, cesm, cmip, erai, and wrf directories will require the mygis packages as well.
 The Python script `ideal_linear.py` will require Nio to be installed with `pip install nio`.
 
-## Install With Conda
+## Setup Environment
+### Install With Conda
 ```bash
 $ conda env create -f environment.yml --prefix /path/to/install/icar_env
 $ conda activate icar_env
+$ conda env config vars set PYTHONPATH=$(pwd)/lib:$PYTHONPATH
+
+reactivate environment, this will be saved for future use
+$ conda activate icar_env
 ```
 
-## Install With Pip
+### Install With Pip
 ```bash
 $ pip install -r requirements.txt
 ```
+Make sure the `lib` directory is in the `PYTHONPATH`, add to `.bashrc` or other startup files for repeat use.
+```bash
+$ export PYTHONPATH=$(pwd)/lib:$PYTHONPATH
+```
+
+
+## Deprecated Scripts
+The units package used in the `cmip/convert.py, cesm/bias_correct.py` and `gen_sounding.py` scripts has been deprecated.

helpers/ccsm/config.py

@@ -15,28 +15,28 @@ def set_bounds(info):
     ccsm_file=atm_file.replace("_Y_","2006").replace("_M_","01").replace("_D_","01").replace("_VAR_","hus")
     ccsm_file=glob.glob(ccsm_file)[0]
     varlist=["lat","lon"]
-    
+
     lat=io.read_nc(ccsm_file,varlist[0]).data
     lon=io.read_nc(ccsm_file,varlist[1]).data-360
-    
+
     info.xmin=np.where(lon>=info.lon[0])[0][0]
     info.xmax=np.where(lon<=info.lon[1])[0][-1]+1
     info.ymin=np.where(lat>=info.lat[0])[0][0]
     info.ymax=np.where(lat<=info.lat[1])[0][-1]+1
-    
+
     lon,lat=np.meshgrid(lon[info.xmin:info.xmax],lat[info.ymin:info.ymax])
     info.lat_data=lat
     info.lon_data=lon
-    
+
 def make_timelist(info,hrs=6.0):
     dt=datetime.timedelta(hrs/24)
-    info.ntimes=np.int(np.round((info.end_date-info.start_date).total_seconds()/60./60./hrs))
+    info.ntimes=np.int64(np.round((info.end_date-info.start_date).total_seconds()/60./60./hrs))
     info.times=[info.start_date+dt*i for i in range(info.ntimes)]
 
 def update_info(info):
     make_timelist(info)
     set_bounds(info)
-    
+
 
 def parse():
     parser= argparse.ArgumentParser(description='Convert CCSM files to ICAR input forcing files')
@@ -51,24 +51,24 @@ def parse():
     parser.add_argument('sfcdir',    nargs="?",action='store',help="CCSM surface data file location",    default="ccsm_sfc/")
     parser.add_argument('atmfile',   nargs="?",action='store',help="CCSM atmospheric files",             default="_VAR__6hrLev_CCSM4_rcp85_r6i1p1__Y__M__D_00-*.nc")
     parser.add_argument('sfcfile',   nargs="?",action='store',help="CCSM surface files",                 default="_VAR__3hr_CCSM4_rcp85_r6i1p1__Y__M__D_0000-*.nc")
-    
+
     parser.add_argument('-v', '--version',action='version',
             version='CCSM2ICAR v'+version)
     parser.add_argument ('--verbose', action='store_true',
             default=False, help='verbose output', dest='verbose')
     args = parser.parse_args()
-    
+
     date0=args.start_date.split("-")
     start_date=datetime.datetime(int(date0[0]),int(date0[1]),int(date0[2]))
 
     date0=args.end_date.split("-")
     end_date=datetime.datetime(int(date0[0]),int(date0[1]),int(date0[2]))
-    
+
     info=Bunch(lat=[float(args.lat_s),float(args.lat_n)],
                lon=[float(args.lon_w),float(args.lon_e)],
                start_date=start_date,         end_date=end_date,
                atmdir=args.dir+args.atmdir,   sfcdir=args.dir+args.sfcdir,
                atmfile=args.atmfile,          sfcfile=args.sfcfile,
                version=version)
-    
+
     return info

helpers/cmip/io_routines.py

@@ -17,9 +17,9 @@ def read_nc(filename,var="data",proj=None,returnNCvar=False):
         data:raw data as an array
         proj:string representation of the projection information
         atts:data attribute dictionary (if any)
-    if (returnNCvar==True) then the netCDF4 file is note closed and the netCDF4 
-        representation of the variable is returned instead of being read into 
-        memory immediately.  
+    if (returnNCvar==True) then the netCDF4 file is note closed and the netCDF4
+        representation of the variable is returned instead of being read into
+        memory immediately.
     '''
     d=netCDF4.Dataset(filename, mode='r',format="nc")
     outputdata=None
@@ -39,8 +39,8 @@ def read_nc(filename,var="data",proj=None,returnNCvar=False):
     if proj!=None:
         projection=d.variables[proj]
         outputproj=str(projection)
-    
-    
+
+
     if returnNCvar:
         return Bunch(data=outputdata,proj=outputproj,ncfile=d,atts=attributes)
     d.close()
@@ -53,7 +53,7 @@ def find_atm_file(time,varname,info):
     file_base= file_base.replace("_Y_",str(time.year))
     file_base= file_base.replace("_EXP_",info.experiment)
     atm_file = file_base.replace("_ENS_",info.ensemble)
-    
+
     print(atm_file)
     filelist = glob.glob(atm_file)
     filelist.sort()
@@ -67,7 +67,7 @@ def find_sst_file(time,info):
     file_base= file_base.replace("_Y_",str(time.year))
     file_base= file_base.replace("_EXP_",info.experiment)
     sst_file = file_base.replace("_ENS_",info.ensemble)
-    
+
     print(sst_file)
     filelist=glob.glob(sst_file)
     filelist.sort()
@@ -76,7 +76,7 @@ def find_sst_file(time,info):
 
 def load_atm(time,info):
     """Load atmospheric variable from a netcdf file"""
-    
+
     outputdata=Bunch()
 
     for s,v in zip(icar_atm_var,atmvarlist):
@@ -98,24 +98,24 @@ def load_atm(time,info):
             outputdata[varname]=np.concatenate([outputdata[varname],newdata])
         else:
             outputdata[varname]=newdata
-        
+
         varname="p"
         newdata=info.read_pressure(atmfile)[:,:,info.ymin:info.ymax,info.xmin:info.xmax]
         if varname in outputdata:
             outputdata[varname]=np.concatenate([outputdata[varname],newdata])
         else:
             outputdata[varname]=newdata
-            
+
         outputdata.ntimes = outputdata.p.shape[0]
-        
+
     # outputdata.times=info.read_time(atmfile)
     try:
         calendar = mygis.read_attr(atmfile_list[0], "calendar", varname="time")
-    except KeyError,IndexError:
+    except (KeyError,IndexError):
         calendar = None
-    
+
     outputdata.calendar = calendar
-    
+
     return outputdata
 
 def load_sfc(time,info):
@@ -137,5 +137,3 @@ def load_data(time,info):
     atm=load_atm(time,info)
     sfc=load_sfc(time,info)
     return Bunch(sfc=sfc,atm=atm)
-
-

helpers/erai/config.py

@@ -3,9 +3,8 @@
 
 import numpy as np
 
-from bunch import Bunch
 import sys
-sys.path.append('../lib') # required to find mygis.py. Update as required.
+from bunch import Bunch
 import mygis
 
 import io_routines as io
@@ -53,7 +52,7 @@ def set_bounds(info):
 def make_timelist(info):
     hrs=6.0
     dt=datetime.timedelta(hrs/24)
-    info.ntimes=np.int(np.round((info.end_date-info.start_date).total_seconds()/60./60./hrs))
+    info.ntimes=np.int64(np.round((info.end_date-info.start_date).total_seconds()/60./60./hrs))
     info.times=[info.start_date+dt*i for i in range(info.ntimes)]
 
 def update_info(info):

helpers/gen_bc.py

@@ -26,8 +26,8 @@ def update_base(base,filename,nz):
 def main():
     filename="bc"
     nx,ny,nz,nt=(20.,20,10,24)
-    dims=[nt,nz,ny,nx]
-    
+    dims=[nt,nz,ny,int(nx)]
+
     lonmin=-110.0; lonmax=-100.0; dlon=(lonmax-lonmin)/nx
     latmin=35.0; latmax=45.0; dlat=(latmax-latmin)/ny
 
@@ -41,32 +41,32 @@ def main():
         update_base(base,"sounding.txt",nz)
         nz=base.th.size
         dims=[nt,nz,ny,nx]
-    
+
     u=np.zeros(dims,dtype="f")+base.u
     w=np.zeros(dims,dtype="f")+base.w
     v=np.zeros(dims,dtype="f")+base.v
     qv=np.zeros(dims,dtype="f")+base.qv
     qc=np.zeros(dims,dtype="f")+base.qc
     coscurve=np.cos(np.arange(dims[2])/dims[2]*2*np.pi+np.pi)+1
-    hgt=(coscurve*1000).reshape((1,nx)).repeat(ny,axis=0)
-    
+    hgt=(coscurve*1000).reshape((1,int(nx))).repeat(ny,axis=0)
+
     lon=np.arange(lonmin,lonmax,dlon)
     lat=np.arange(latmin,latmax,dlat)
     lon,lat=np.meshgrid(lon,lat)
-    
+
     dz=np.zeros(dims)+base.dz
-    z=np.zeros(dims,dtype="f")+base.z.reshape((1,nz,1,1))+hgt.reshape((1,1,ny,nx))
-    
+    z=np.zeros(dims,dtype="f")+base.z.reshape((1,nz,1,1))+hgt.reshape((1,1,ny,int(nx)))
+
     layer1=(dz[0,0,:,:]/2)
     z[0,0,:,:]+=layer1
     for i in range(1,int(nz)):
         z[:,i,:,:]=z[:,i-1,:,:]+(dz[:,i-1,:,:]+dz[:,i,:,:])/2.0
-    
+
     p=np.zeros(dims,dtype="f")+base.p
     adjust_p(p,0.0,z)
     th=np.zeros(dims,dtype="f")+base.th
-    
-    
+
+
     d4dname=("t","z","y","x")
     d3dname=("z","y","x")
     d2dname=("y","x")
@@ -86,7 +86,7 @@ def main():
     if fileexists:
         print("Removing : "+fileexists[0])
         os.remove(fileexists[0])
-    
+
     io.write(filename,  u,varname="U", dims=d4dname,dtype="f",attributes=dict(units="m/s",description="Horizontal (x) wind speed"),
             extravars=othervars)
 

helpers/gen_ideal.py

@@ -30,7 +30,7 @@ def adjust_p(p,h,dz):
 
 def update_base(base,filename,nz):
     """update the base information using data from a sounding file
-    
+
     filename should be a space delimited text file with 3 columns
         height [m], potential temperature [K], and specific humidity [g/kg]"""
     print("Using Sounding from : "+filename)
@@ -43,7 +43,7 @@ def update_base(base,filename,nz):
 
 def build_topography(experiment,dims):
     """create the topography to be used for a given case study"""
-    
+
     # experiment                        D1      D2      D3
     # h         height of the hill    1800    1400    1040   [meters]
     # sigma     half-width              60      40       3.1 [grid cells]
@@ -60,8 +60,8 @@ def build_topography(experiment,dims):
     x     = np.linspace(0,Lx,Nx)                  # % distance array (m)
     zo    = [1700.0,2000.0,2200.0][experiment]    # mountain base height (m) NOT REALLY USED CORRECTLY YET, NEED to truncate the sounding...
     zo    = 0.0
-    
-    
+
+
     zs=hm/(1.0+((x/dx-xm)/am)**2.)
     # zs-=zs[0]
     zs=zs.reshape((1,dims[3])).repeat(dims[2],axis=0)
@@ -71,8 +71,8 @@ def main():
     filename="ideal_{}_{}".format(case_study,int(wind_speed))
     print(filename)
     nx,nz,ny=master_dims[case_study]
-    dims=[1,nz,ny,nx]
-    
+    dims=[1,int(nz),int(ny),int(nx)]
+
     # this is just arbitrary for now
     dlon=dx/111.1
     dlat=dx/111.1
@@ -89,7 +89,7 @@ def main():
         update_base(base,"sounding.txt",nz)
         nz=base.th.size
         dims=[1,nz,ny,nx]
-    
+
     udims=copy(dims)
     udims[-1]+=1
     vdims=copy(dims)
@@ -99,40 +99,40 @@ def main():
     v=np.zeros(vdims,dtype="f")+base.v
     qv=np.zeros(dims,dtype="f")+base.qv
     qc=np.zeros(dims,dtype="f")+base.qc
-    
+
     # simple topography = a cosine
     # coscurve=np.cos(np.arange(dims[3])/dims[3]*2*np.pi+np.pi)+1
     # hgt=(coscurve*1000).reshape((1,nx)).repeat(ny,axis=0)
     hgt=build_topography(case_study,dims)
-    
-    lon=np.arange(lonmin,lonmax,dlon)[:nx]
+
+    lon=np.arange(lonmin,lonmax,dlon)[:int(nx)]
     lat=np.arange(latmin,latmax,dlat)[:ny]
     lon,lat=np.meshgrid(lon,lat)
 
-    ulon=np.arange(lonmin-dlon/2,lonmax+dlon/2,dlon)[:nx+1]
+    ulon=np.arange(lonmin-dlon/2,lonmax+dlon/2,dlon)[:int(nx)+1]
     ulat=np.arange(latmin,latmax,dlat)[:ny]
     ulon,ulat=np.meshgrid(ulon,ulat)
 
-    vlon=np.arange(lonmin,lonmax,dlon)[:nx]
+    vlon=np.arange(lonmin,lonmax,dlon)[:int(nx)]
     vlat=np.arange(latmin-dlat/2,latmax+dlat/2,dlat)[:ny+1]
     vlon,vlat=np.meshgrid(vlon,vlat)
-    
+
     dz=np.zeros(dims)+base.dz
-    z=np.zeros(dims,dtype="f")+base.z.reshape((1,nz,1,1))+hgt.reshape((1,1,ny,nx))
-    
+    z=np.zeros(dims,dtype="f")+base.z.reshape((1,nz,1,1))+hgt.reshape((1,1,ny,int(nx)))
+
     layer1=(dz[0,:,:]/2)
     z[0,:,:]+=layer1
     for i in range(1,int(nz)):
         z[:,i,:,:]=z[:,i-1,:,:]+(dz[:,i-1,:,:]+dz[:,i,:,:])/2.0
-    
+
     p=np.zeros(dims,dtype="f")+base.p
     adjust_p(p,0.0,z)
     th=np.zeros(dims,dtype="f")+base.th
-    
-    lat=lat.reshape((1,ny,nx))
-    lon=lon.reshape((1,ny,nx))
-    hgt=hgt.reshape((1,ny,nx))
-    
+
+    lat=lat.reshape((1,ny,int(nx)))
+    lon=lon.reshape((1,ny,int(nx)))
+    hgt=hgt.reshape((1,ny,int(nx)))
+
     d3dname=("t","z","y","x")
     ud3dname=("t","z","y","xu")
     ud2dname=("t","y","xu")
@@ -164,13 +164,12 @@ def main():
     if fileexists:
         print("Removing : "+fileexists[0])
         os.remove(fileexists[0])
-    
+
     io.write(filename,  u,varname="U", dims=ud3dname,dtype="f",attributes=dict(units="m/s",description="Horizontal (x) wind speed"),
             extravars=othervars)
 
 
 if __name__ == '__main__':
-    global wind_speed, case_study
     for case in range(3):
         for ws in [5,10,15,25]:
             wind_speed=float(ws)