Changeset 1974 for TI02-CSML/trunk


Ignore:
Timestamp:
09/01/07 10:48:16 (13 years ago)
Author:
domlowe
Message:

improved determineCRS method

File:
1 edited

Legend:

Unmodified
Added
Removed
  • TI02-CSML/trunk/csml/csmllibs/csmlcrs.py

    r1972 r1974  
    3737        self.systems['ndg:crs:xyt']=crs 
    3838         
     39        #define unknown 3D CRS: 
     40        crs=CRSystem(srsName='ndg:crs:unknown3d', axes=['unknown','unknown','unknown']) 
     41        self.systems['ndg:crs:unknown3d']=crs 
     42         
     43        #define unknown 4D CRS: 
     44        crs=CRSystem(srsName='ndg:crs:unknown4d', axes=['unknown','unknown','unknown','unknown']) 
     45        self.systems['ndg:crs:unknown4d']=crs 
     46         
    3947        #define unknown CRS: 
    4048        crs=CRSystem(srsName='ndg:crs:unknown', axes=['unknown']) 
     
    4957     
    5058    def determineCRS(self, axes, units): 
    51         #given any list of axis names and a list of units for these axes attempt to determine the CRS and return the CRSystem object 
    52         crs=self.systems['ndg:crs:unknown'] 
     59        '''given any list of axis names and a list of units for these axes attempt to determine the CRS and return the CRSystem object and an axis order list. 
     60        e.g passing in: (axes=['t', 'ht', 'latitude', 'longitude'],units=['days since 1991-09-01 00:00:00', 'm', 'degrees_north', 'degrees_east']) 
     61        returns self.systems['ndg:crs:xypt'], [3,2,1,0] 
     62        Assumes the units are ordered to correspond with the axes! 
     63        ''' 
     64        axisorder=[] 
    5365        if len(axes)==3: 
    54             #it's a 3d crs 
    55             for axis in axes: 
    56                 for unit in units: 
    57                     if string.lower(unit) in ['second', 'seconds', 's', 'mins','minute','minutes','hour','hours','h','hr','hrs','day','days']: 
    58                         #unit is time 
    59                         return self.systems['ndg:crs:xyt'] 
     66            crs=self.systems['ndg:crs:unknown3d'] 
    6067        elif len(axes)==4: 
    61             #it's a 4d crs 
    62             pass        
    63         return crs 
     68            crs=self.systems['ndg:crs:unknown4d'] 
     69        else : 
     70            crs=self.systems['ndg:crs:unknown'] 
     71        crsMap=[] 
     72        #this can be extended to accomodate more units 
     73        for axis in axes: 
     74            print axis 
     75            unit=units[axes.index(axis)] 
     76            print unit 
     77            unittype='unknown' 
     78            if string.lower(unit) in ['second', 'seconds', 's', 'mins','minute','minutes','hour','hours','h','hr','hrs','day','days']: 
     79                unittype='Time' 
     80            elif string.lower(unit)[:10] in ['days since', 'seconds si', 'minutes si', 'hours sinc']: 
     81                unittype='Time' 
     82            elif string.lower(unit) in ['mbar', 'pa','level']: 
     83                unittype='Pressure' 
     84            elif string.lower(unit) in ['m', 'km']: 
     85                unittype='Height' 
     86            elif string.lower(unit) in ['degrees_north']: 
     87                unittype='Lat' 
     88            elif string.lower(unit) in ['degrees_east']: 
     89                unittype='Lon' 
     90            crsMap.append(unittype) 
     91        print 'CRS_MAP %s'%crsMap 
     92         
     93        #now try and match up crsMap with known crsystems.         
     94        for system in self.systems: 
     95            if len(self.systems[system].axes) != len(crsMap): 
     96                #wrong dimensionality, skip this crs 
     97                continue 
     98            match=0 
     99            axisorder=[] 
     100            for ax in self.systems[system].axes: 
     101                if ax in crsMap: 
     102                    axisorder.append(crsMap.index(ax)) 
     103                else: 
     104                    match = 1 
     105                    break 
     106            if match == 0:  
     107                crs=self.systems[system] 
     108                break 
     109                         
     110        return crs, axisorder 
    64111         
    65112def main(): 
     
    73120    print 'AXES: %s'%axs 
    74121    print 'UNITS: %s'%uns 
    75     crs= cat.determineCRS(axes=axs, units=uns) 
     122    crs, axisorder= cat.determineCRS(axes=axs, units=uns) 
    76123    print 'srsName = %s' %crs.srsName 
    77124    print 'srsDimension = %s' %crs.srsDimension 
     
    85132    print 'AXES: %s'%axs 
    86133    print 'UNITS: %s'%uns 
    87     crs=cat.determineCRS(axes=axs, units=uns) 
     134    crs, axisorder=cat.determineCRS(axes=axs, units=uns) 
    88135    print 'srsName = %s' %crs.srsName 
    89136    print 'srsDimension = %s' %crs.srsDimension 
     
    97144        v=f.variables[var] 
    98145        axes=v.getAxisIds() 
    99  #   print axes 
    100146        units =[] 
    101147        for axis in axes: 
    102148            ax=f.getAxis(axis) 
    103149            units.append(ax.attributes['units']) 
    104 #    print units 
    105150        crsdictionary[str(axes)]=units 
    106     #print crsdictionary 
    107151    n=2 
    108152    #test each definition in dictionary 
     
    111155        print 'TEST %d'%n 
    112156        print '*************************************' 
    113         axs =key 
     157       #make list from string. 
     158        axs=[] 
     159        for item in key.split(): 
     160            axs.append(item) 
    114161        uns=crsdictionary[key] 
    115162        print 'AXES: %s'%axs 
    116163        print 'UNITS: %s'%uns 
    117         crs= cat.determineCRS(axes=axs, units=uns) 
     164        crs, axisorder= cat.determineCRS(axes=axs, units=uns) 
    118165        print 'srsName = %s' %crs.srsName 
    119166        print 'srsDimension = %s' %crs.srsDimension 
    120167        print 'axisLabels =%s' %crs.axisLabels 
     168        print 'axis order =%s'%axisorder 
    121169        print '**************************************' 
    122170   #test getting CRS by name: 
Note: See TracChangeset for help on using the changeset viewer.