Changeset 690 in lmdz_wrf for trunk/tools/nc_var_tools.py


Ignore:
Timestamp:
Feb 4, 2016, 5:48:14 PM (9 years ago)
Author:
lfita
Message:

Changing the `check_arguments' functions

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/tools/nc_var_tools.py

    r689 r690  
    126126    return vals
    127127
    128 def check_arguments(funcname,Nargs,args,char,expectargs):
     128def check_arguments(funcname,args,expectargs,char):
    129129    """ Function to check the number of arguments if they are coincident
    130130    check_arguments(funcname,Nargs,args,char)
    131131      funcname= name of the function/program to check
    132       Nargs= theoretical number of arguments
    133132      args= passed arguments
     133      expectargs= expected arguments
    134134      char= character used to split the arguments
    135135    """
    … …  
    138138
    139139    Nvals = len(args.split(char))
     140    Nargs = len(expectargs.split(char))
     141
    140142    if Nvals != Nargs:
    141143        print errormsg
    142144        print '  ' + fname + ': wrong number of arguments:',Nvals," passed to  '",   \
    143145          funcname, "' which requires:",Nargs,'!!'
    144         print '    given arguments:',args.split(char)
    145         print '    expected arguments:',expectargs
     146        print '     # given expected _______'
     147        Nmax = np.max([Nvals, Nargs])
     148        for ia in range(Nmax):
     149            if ia > Nvals-1:
     150                aval = ' '
     151            else:
     152                aval = args.split(char)[ia]
     153            if ia > Nargs-1:
     154                aexp = ' '
     155            else:
     156                aexp = expectargs.split(char)[ia]
     157
     158            print '    ', ia, aval, aexp
    146159        quit(-1)
    147160
    … …  
    76467659      'newlatdim]|[oldlonname,oldlatname]|[oldlondim,oldlatdim]'
    76477660 
    7648     check_arguments(fname,len(expectargs.split('|')),values,'|',expectargs)
     7661    check_arguments(fname,values,expectargs,'|')
    76497662
    76507663    newprojectionfile = values.split('|')[0]
    … …  
    1121511228    arguments = '[trajfile]@[Tbeg],[lonname],[latname],[zname],[timename],' +        \
    1121611229      '[timekind],[boxsize],[circler]]'
    11217     check_arguments(fname,len(arguments.split(',')),arguments,',',8)
     11230    check_arguments(fname,values,arguments,',')
    1121811231
    1121911232    trajfile = values.split(',')[0].split('@')[0]
    … …  
    1368013693
    1368113694    expectargs = '[sectionname],[kout]'
    13682     check_arguments(fname,len(expectargs.split(',')),values,',',expectargs)
     13695    check_arguments(fname,values,expectargs,',')
    1368313696
    1368413697    secname = values.split(',')[0]
    … …  
    1385513868    expectargs = expectargs + '[geogres]'
    1385613869 
    13857     check_arguments(fname,len(expectargs.split(',')),values,',',expectargs)
     13870    check_arguments(fname,values,expectargs,',')
    1385813871
    1385913872    namelist = values.split(',')[0]
    … …  
    1399014003    expectargs = '[Ngrid],[dimnh1],[dimnh2],[dimvnh1],[dimvnh2]'
    1399114004 
    13992     check_arguments(fname,len(expectargs.split(',')),values,',',expectargs)
     14005    check_arguments(fname,values,expectargs,',')
    1399314006
    1399414007    Ngrid = int(values.split(',')[0])
    … …  
    1417714190    expectargs = '[kind],[kindn],[netcdffileref]'
    1417814191 
    14179     check_arguments(fname,len(expectargs.split(',')),values,',',expectargs)
     14192    check_arguments(fname,values,expectargs,',')
    1418014193
    1418114194    kind = values.split(',')[0]
    … …  
    1430414317    expectargs = '[tdimn],[tvarn]'
    1430514318 
    14306     check_arguments(fname,len(expectargs.split(',')),values,',',expectargs)
     14319    check_arguments(fname,values,expectargs,',')
    1430714320
    1430814321    tdimn = values.split(',')[0]
    … …  
    1440614419    expectargs = '[dimensions]|[varattributes]|[kind]'
    1440714420 
    14408     check_arguments(fname,len(expectargs.split('|')),values,'|',expectargs)
     14421    check_arguments(fname,values,expectargs,'|')
    1440914422
    1441014423    dimensions = values.split('|')[0].split(',')
    … …  
    1449414507    expectargs = '[dimensions]|[varattributes]|[kind]'
    1449514508 
    14496     check_arguments(fname,len(expectargs.split('|')),values,'|',expectargs)
     14509    check_arguments(fname,values,expectargs,'|')
    1449714510
    1449814511    dimensions = values.split('|')[0].split(',')
    … …  
    1465114664    expectargs = '[dimname]:[size]:[position]:[unlimited]'
    1465214665 
    14653     check_arguments(fname,len(expectargs.split(':')),values,':',expectargs)
     14666    check_arguments(fname,values,expectargs,':')
    1465414667
    1465514668    dimname = values.split(':')[0]
    … …  
    1546215475    arguments = '[trajfile]@[Tbeg],[lonname],[latname],[zname],[timename],' +        \
    1546315476      '[timekind],[Nangle],[radii]]'
    15464     check_arguments(fname,len(arguments.split(',')),arguments,',',8)
     15477    check_arguments(fname,values,arguments,',')
    1546515478
    1546615479    trajfile = values.split(',')[0].split('@')[0]
    … …  
    1629016303
    1629116304    arguments = '[dates],[fold]'
    16292     check_arguments(fname, len(values.split(',')), arguments, ',',                   \
    16293       len(arguments.split(',')))
     16305    check_arguments(fname, values, arguments, ',')
    1629416306
    1629516307    dates = values.split(',')[0].split(':')
    … …  
    1638416396
    1638516397    arguments = '[tdim]:[tvdim]:[seltimes]'
    16386     check_arguments(fname, len(values.split(':')), arguments, ':',                   \
    16387       len(arguments.split(',')))
     16398    check_arguments(fname, values, arguments, ':')
    1638816399
    1638916400    tdim = values.split(':')[0]
    … …  
    1678516796
    1678616797    arguments = '[longitude]:[latitude]:[dnames]:[vdnames]'
    16787     check_arguments(fname, len(values.split(':')), arguments, ':',                   \
    16788       len(arguments.split(',')))
     16798    check_arguments(fname, values, arguments, ':')
    1678916799
    1679016800    longitude = np.float(values.split(':')[0])
    … …  
    1728917299    if weightk == 'varnames':
    1729017300        arguments = '[weightk],[varname],[oper],[xdimname],[ydimname],[addvars]'
    17291         check_arguments(fname, len(values.split(',')), arguments, ',',               \
    17292           len(arguments.split(',')))
     17301        check_arguments(fname, values, arguments, ',')
    1729317302        varname = values.split(',')[1]
    1729417303        oper = values.split(',')[2]     
    … …  
    1729917308            quit(-1)
    1730017309    elif weightk == 'reglonlat':
    17301         arguments = '[weightk],[lonname],[latname],[xdimname],[ydimname]'
    17302         check_arguments(fname, len(values.split(',')), arguments, ',',               \
    17303           len(arguments.split(',')))
     17310        arguments = '[weightk],[lonname],[latname],[xdimname],[ydimname],[addvars]'
     17311        check_arguments(fname, values, arguments, ',')
    1730417312        lonname = values.split(',')[1]
    1730517313        latname = values.split(',')[2]
    … …  
    1748817496
    1748917497        ivarwgtv = np.abs(np.cos(latv*np.pi/180.))/(lonv.shape[0]*lonv.shape[1])
    17490         TOTsumwgt = np.sum(ivarwgt)
     17498        TOTsumwgt = np.sum(ivarwgtv)
    1749117499        if len(loopshape) == 1:
    1749217500            newvals = np.ones((loopshape[0]), dtype=np.float)*fillValueF
    … …  
    1773417742    if project == 'lonlat':
    1773517743        Nlon = int(2. * np.pi * EarthR / dxres) + 1
    17736         Nlat = int(np.pi * EarthR / dxres) + 1
     17744        Nlat = int(np.pi * EarthR / dyres) + 1
    1773717745
    1773817746        lonProj = np.zeros((Nlon), dtype=np.float)
    … …  
    1775017758        lonProj = []
    1775117759        latProj = []
     17760        lonProjcirc = []
     17761        latProjcirc = []
     17762        Nlonmaxcirc = -10.
     17763
    1775217764        lonProj.append(0.)
    1775317765        latProj.append(-90.)
     17766        lonProjcirc.append(0.)
     17767        latProjcirc.append(-90.)
    1775417768        # Differential on latitudes
    17755         latarc = dyres / (np.pi * EarthR)
     17769        Nlat = int(np.pi * EarthR / dyres) + 3
     17770        latarc = np.pi / (Nlat - 1)
     17771
     17772        Nloncirc = np.zeros((Nlat), dtype=int)
     17773        Nloncirc[0] = 1
    1775617774
    1775717775        lat=-np.pi/2.
     17776        ilat = 1
    1775817777        while lat + latarc <= np.pi/2.:
    1775917778            lon=0.
    1776017779            lat = lat + latarc
    1776117780            lonarc = dxres / (EarthR * np.cos(lat))
     17781            lonscirc = []
     17782            latProjcirc.append(lat*RadDeg)
    1776217783            while lon <= 2.*np.pi:
    1776317784#                print lonarc, latarc, lon, lat, [lon*RadDeg, lat*RadDeg]
    17764                 lonProj.append(-180. + lon*RadDeg) 
     17785                lonProj.append(-180. + lon*RadDeg)
    1776517786                latProj.append(lat*RadDeg)
     17787                lonscirc.append(-180. + lon*RadDeg)
    1776617788                lon = lon + lonarc
    1776717789                Npts = Npts + 1
     17790            lonProjcirc.append(lonscirc)
     17791            Nloncirc[ilat] = len(lonscirc)
     17792
     17793            if Nloncirc[ilat] > Nlonmaxcirc: Nlonmaxcirc = Nloncirc[ilat]
     17794            ilat = ilat + 1
     17795
    1776817796        lonProj.append(0.)
    1776917797        latProj.append(90.)
     17798        lonProjcirc.append(0.)
     17799        latProjcirc.append(90.)
     17800        Nloncirc[ilat] = 1
    1777017801        Npts = Npts + 1
    1777117802
    … …  
    1778817819        newdim = onewnc.createDimension('Npts', Npts)
    1778917820        newdim = onewnc.createDimension('Ncoord', 2)
     17821        newdim = onewnc.createDimension('xcirc', Nlonmaxcirc)
     17822        newdim = onewnc.createDimension('ycirc', Nlat)
    1779017823
    1779117824# Variables
    1779217825        newvar = onewnc.createVariable('lon', 'f8', ('x'))
    17793         basicvardef(newvar, 'longitude', 'Longitude', 'degress_East')
     17826        basicvardef(newvar, 'longitude', 'Longitude', 'degrees_East')
    1779417827        newvar[:] = lonProj
    1779517828        newvar.setncattr('axis', 'X')
    1779617829        newvar.setncattr('_CoordinateAxisType', 'Lon')   
    1779717830
     17831# Longitudes sorted by latitudes
     17832        if len(Npoints.shape) == 1:
     17833            newvar = onewnc.createVariable('Nloncirc', 'i', ('ycirc'))
     17834            newvar[:] = Nloncirc
     17835            basicvardef(newvar, 'Nloncirc', 'Number of longitudes by latitude circs', '-')
     17836
     17837            newvar = onewnc.createVariable('loncirc', 'f8', ('ycirc','xcirc'))
     17838            basicvardef(newvar, 'longitude', 'Longitude by latitude circs', 'degrees_East')
     17839            for iL in range(Nlat):
     17840                newvar[iL,0:Nloncirc[iL]] = np.asarray(lonProjcirc[iL])
     17841
    1779817842        newvar = onewnc.createVariable('lat', 'f8', ('y'))
    17799         basicvardef(newvar, 'latitude', 'Latitude', 'degress_North')
     17843        basicvardef(newvar, 'latitude', 'Latitude', 'degrees_North')
    1780017844        newvar[:] = latProj
    1780117845        newvar.setncattr('axis', 'Y')
    1780217846        newvar.setncattr('_CoordinateAxisType', 'Lat')   
     17847
     17848        newvar = onewnc.createVariable('latcirc', 'f8', ('ycirc'))
     17849        basicvardef(newvar, 'latitude', 'Latitude', 'degrees_North')
     17850        newvar[:] = latProjcirc
    1780317851
    1780417852        if len(Npoints.shape) == 2:
    … …  
    1781017858        newvar[:] = Npoints
    1781117859        newvar.setncattr('coordinates',  'lon lat')
    17812         print Npoints.shape
    17813    
    1781417860
    1781517861# Global attributes
    … …  
    1782717873        print fname + "Successfull written of '" + ofile + "' !!'"
    1782817874
    17829     return lonProj, latProj
    17830 
    17831 #lonlatProj(100000.,100000.,'lonlat',True)
     17875    if len(Npoints.shape) == 2:
     17876        return lonProj, latProj
     17877    else:
     17878        return lonProjcirc, latProj
     17879
     17880#lonlatProj(100000.,100000.,'lonlat_dxdyFix',True)
    1783217881#quit()
    1783317882
    … …  
    1784817897    """
    1784917898    import numpy.ma as ma
     17899    import subprocess as sub
     17900
    1785017901    fname = 'Partialmap_Entiremap'
    1785117902
    … …  
    1785617907
    1785717908    ofile = 'EntireGlobalMap.nc'
     17909
     17910    lonlatProjfile = filen.split('.')[0] + '_lonlatProj.nc'
    1785817911
    1785917912    onc = NetCDFFile(filen, 'r')
    … …  
    1786717920    latname = values.split(',')[1]
    1786817921    fval0 = values.split(',')[2]
    17869     resoltuion = np.float(values.split(',')[3])
     17922    resolution = np.float(values.split(',')[3])
    1787017923
    1787117924    if not onc.variables.has_key(lonname):
    … …  
    1789217945    minlon = np.min(np.abs(lonvs))
    1789317946
     17947#    mindiff = np.min([minlat, minlon])/1.
     17948    mindiff = np.min([minlat, minlon])*10.e6
     17949
    1789417950    print '  ' + fname + ': Closest latitude to Equator:', minlat
    1789517951    print '  ' + fname + ': Closest longitude to Greenwich Meridian:', minlon
     17952    print '  ' + fname + ': Minium distance for point coincidence:', mindiff
    1789617953
    1789717954# Existing longitudes along/closest to the Equator
    … …  
    1792617983        mindgrlat = None
    1792717984
    17928 # Global sorting
    17929 
    17930 #    TOTlons = lonvs.shape[0]
    17931 #    TOTlats = lonvs.shape[0]
    17932 
    17933 #    sortlons = np.sort(lonvs)
    17934 #    sortlats = np.sort(latvs)
    17935 
    17936 #    difflons = sortlons[1:TOTlons-1] - sortlons[0:TOTlons-2]
    17937 #    difflats = sortlats[1:TOTlons-1] - sortlats[0:TOTlons-2]
    17938 
    17939 #    madifflons = ma.masked_equal(difflons, 0.)
    17940 #    madifflats = ma.masked_equal(difflats, 0.)
    17941 
    17942 #    print 'Nmadifflons:',len(madifflons.compressed()),'Nmadifflats:',len(madifflats.compressed())
    17943 
    17944 #    sortdifflons = np.sort(madifflons)
    17945 #    sortdifflats = np.sort(madifflats)
    17946 
    17947 #    print 'sortdifflons:',sortdifflons[0:10]
    17948 #    print 'sortdifflats:',sortdifflats[0:10]
    17949 
    1795017985# Fixing in case it has not worked
    1795117986    if mindeqlon is not None and mindgrlat is None: mindgrlat = mindeqlon
    … …  
    1795617991        quit(-1)
    1795717992
    17958     if nlon < 0.:
    17959         lonmap = np.arange(-180.+mindeqlon,180.,mindeqlon)
     17993    if not os.path.isfile(lonlatProjfile):
     17994        print lonlatProjfile
     17995        quit(-1)
     17996        lonmap, latmap = lonlatProj(resolution, resolution, 'lonlat_dxdyFix', True)
     17997        sub.call(['mv','lonlatProj.nc',lonlatProjfile], shell=True)
     17998        Ntotvals = len(lonmap)
    1796017999    else:
    17961         lonmap = np.arange(0.,360.-mindeqlon,mindeqlon)
    17962 
    17963     if nlat < 0.:
    17964         latmap = np.arange(-90.+mindgrlat,90.,mindgrlat)
    17965     else:
    17966         latmap = np.arange(0.,180.-mindgrlat,mindgrlat)
    17967 
    17968     print '  ' + fname + '_________'
    17969     print '     final longitudes:',lonmap
    17970     print '     final latitudes:',latmap
    17971 
    17972     Nmaplon = len(lonmap)
    17973     Nmaplat = len(latmap)
    17974     print '  ' + fname + ': desired map will be:',Nmaplon,'x',Nmaplat
    17975 
    17976     quit()
     18000        oproj = NetCDFFile(lonlatProjfile, 'r')
     18001        oprojlon = oproj.variables['loncirc']
     18002        oprojlat = oproj.variables['latcirc']
     18003        oNprojlon = oproj.variables['Nloncirc']
     18004        projlat = oprojlat[:]
     18005        Nprojlon = oNprojlon[:]
     18006        Ntotvals = len(oproj.dimensions['Npts'])
     18007        olon = oproj.variables['lon']
     18008        olat = oproj.variables['lat']
     18009        lonmap = olon[:]
     18010        latmap = olat[:]
     18011        omapvals = oproj.variables['points']
     18012
     18013#    Ntotvals = len(lonmap)
     18014    Nmaplon = Ntotvals
     18015    Nmaplat = Ntotvals
    1797718016
    1797818017# Matrix map creation
    … …  
    1798518024    fval = fill_value(vartype, fval0)
    1798618025
    17987     mapval = np.ones((Nmaplat,Nmaplon), dtype=np.float)*fval
    17988 
    17989     Ntotvals = lonvs.shape[0]
     18026# Final map creation
     18027##
     18028    newnc = NetCDFFile(ofile, 'w')
     18029
     18030# dimensions
     18031    newdim = newnc.createDimension('lon',Nmaplon)
     18032    newdim = newnc.createDimension('lat',Nmaplat)
     18033    newdim = newnc.createDimension('Npts', Ntotvals)
     18034
     18035# Variables
     18036    newvar = newnc.createVariable('lon','f8',('lon'))
     18037    basicvardef(newvar, 'lon', 'Longitudes','degrees East')
     18038    newvar[:] = lonmap
     18039    newvar.setncattr('axis', 'X')
     18040    newvar.setncattr('_CoordinateAxisType', 'Lon')   
     18041
     18042    newvar = newnc.createVariable('lat','f8',('lat'))
     18043    basicvardef(newvar, 'lat', 'Latitudes','degrees North')
     18044    newvar[:] = latmap
     18045    newvar.setncattr('axis', 'Y')
     18046    newvar.setncattr('_CoordinateAxisType', 'Lat')   
     18047
     18048# map variable
     18049    Lmapvalsshape = len(omapvals.shape)
     18050    if Lmapvalsshape == 1:
     18051        newvar = newnc.createVariable(varn,vartype,('Npts'), fill_value=fval)   
     18052    else:
     18053        newvar = newnc.createVariable(varn,vartype,('lat','lon'), fill_value=fval)
     18054
     18055    basicvardef(newvar, varn, varlongname, varunits)
     18056    newvar.setncattr('coordinates', 'lon lat')
     18057
     18058    print '  ' + fname + ': re-locating:',Ntotvals,'points...'
    1799018059    for iv in range(Ntotvals):
    17991         if np.mod(iv*100./Ntotvals, 10.) == 0: print '    done:',iv*100./Ntotvals,'%'
    17992         difflon = np.abs(lonmap - lonvs[iv])
    17993         difflat = np.abs(latmap - latvs[iv])
    17994 
    17995 #        print np.min(difflon), ':', np.min(difflat)
    17996 
    17997         ilon = index_mat(difflon, np.min(difflon))
    17998         ilat = index_mat(difflat, np.min(difflat))
    17999         if np.min(difflon) != 0. or np.min(difflat) !=0.:
    18000             print 'No zero diffs!!:', lonvs[iv],latvs[iv],':',np.min(difflon),np.min(difflat),'i j',ilon,ilat
    18001             print lonmap[ilon], latmap[ilat]
     18060#    for iv in range(15):
     18061#        print 'Lluis:',iv,'lon lat:',lonvs[iv],',',latvs[iv]
     18062        difflat = np.abs(projlat - latvs[iv])
     18063        mindiffL = np.min(difflat)
     18064        ilat = index_mat(difflat, mindiffL)
     18065#        print '  Lluis mindiffL:',mindiffL,'<> ilat:',ilat,'lat_ilat:',projlat[ilat],'Nprojlon:',Nprojlon[ilat],'shape oporjlon:',oprojlon[ilat,:].shape, type(oprojlon[ilat,:])
     18066        loncirc = np.zeros((Nprojlon[ilat]), dtype=np.float)
     18067        loncirc[:] = np.asarray(oprojlon[ilat,0:int(Nprojlon[ilat])].flatten())
     18068        difflon = np.abs(loncirc - lonvs[iv])
     18069        mindiffl = difflon.min()
     18070        ilon = index_mat(difflon, mindiffl)
     18071#        print '  difflon:',type(difflon),'shape difflon',difflon.shape,'shape loncirc:', loncirc.shape
     18072#        print '  Lluis mindiffl:',mindiffl,'<> ilon:',ilon,'oprojlon:', loncirc[ilon]
     18073
     18074        percendone(iv,Ntotvals,0.5,'done:')
     18075        if mindiffl > mindiff or mindiffL > mindiff:
     18076            print errormsg
     18077            print '  ' + fname + ': for point #', iv,'lon,lat in incomplet map:',   \
     18078              lonvs[iv], ',', latvs[iv], 'there is not a set of lon,lat in the ' +  \
     18079              'completed map closer than: ', mindiff, '!!'
     18080            print '    minimum difflon:', np.min(difflon), 'difflat:', np.min(difflat)
    1800218081            quit()
    1800318082
    1800418083        if ilon > 0 and ilat > 0:
    18005             mapval[ilat,ilon] = ovar[iv]
     18084#            print 'Lluis shapes:',newvar.shape,'i lon/lat:',ilon,ilat
     18085            if Lmapvalsshape ==1:
     18086                newvar[iv] = ovar[iv]
     18087            else:
     18088                newvar[ilat,ilon] = ovar[iv]
     18089           
     18090            newnc.sync()
    1800618091        else:
    1800718092            print errormsg
    1800818093            print '  ' + fname + ': point:',lonvs[iv],',',latvs[iv],' not relocated !!'
    1800918094            quit(-1)
    18010 
    18011 # Final map creation
    18012 ##
    18013     newnc = NetCDFFile(ofile, 'w')
    18014 
    18015 # dimensions
    18016     newdim = newnc.createDimension('lon',Nmaplon)
    18017     newdim = newnc.createDimension('lat',Nmaplat)
    18018 
    18019 # Variables
    18020     newvar = newnc.createVariable('lon','f8',('lon'))
    18021     basicvardef(newvar, 'lon', 'Longitudes','degrees East')
    18022     nevar[:] = lonmap
    18023 
    18024     newvar = newnc.createVariable('lat','f8',('lat'))
    18025     basicvardef(newvar, 'lat', 'Latitudes','degrees North')
    18026     nevar[:] = latmap
    18027 
    18028 # map variable
    18029     newvar = newnc.createVariable(varn,'f4',('lat','lon'), fill_value=fval)
    18030     basicvardef(newvar, varn, 'Longitudes','degrees East')
    18031     nevar[:] = lonmap
    1803218095
    1803318096    newnc.sync()
    … …  
    1805218115    return
    1805318116
    18054 #Partialmap_Entiremap('longitude,latitude,std', '/home/lluis/etudes/DYNAMICO/ORCHIDEE/interpolation/data/carteveg5km.nc', 'vegetation_map')
     18117#Partialmap_Entiremap('longitude,latitude,std,5000.', '/home/lluis/etudes/DYNAMICO/ORCHIDEE/interpolation/data/carteveg5km.nc', 'vegetation_map')
    1805518118#quit()
    1805618119
Note: See TracChangeset for help on using the changeset viewer.