Changeset 2312 for TI02-CSML/trunk/csml


Ignore:
Timestamp:
22/03/07 15:38:29 (13 years ago)
Author:
domlowe
Message:

adding code to csmlcrs.py for handling unidentified coordinate axis systems, (not fully working yet)

File:
1 edited

Legend:

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

    r2198 r2312  
    1919        # dictionary to hold known CRSystems: 
    2020        self.systems={} 
     21        self.pseudosystems={} 
    2122         
    2223        # define lon lat pressure time CRS:         
     
    5657        crs.timeAxis=1 
    5758         
     59        #define unknown 3D CRS: 
     60        crs=CRSystem(srsName='ndg:crs:undefined2d', axes=['tba','tba']) 
     61        self.pseudosystems['ndg:crs:undefined2d']=crs 
     62         
     63        #define unknown 3D CRS: 
     64        crs=CRSystem(srsName='ndg:crs:undefinedll3d', axes=['longitude','latitude','tba']) 
     65        self.pseudosystems['ndg:crs:undefinedll3d']=crs 
     66        crs.lonAxis=0 
     67        crs.latAxis=1 
     68         
     69        #define unknown 4D CRS: 
     70        crs=CRSystem(srsName='ndg:crs:undefinedllt4d', axes=['longitude','latitude','tba','time']) 
     71        self.pseudosystems['ndg:crs:undefinedllt4d']=crs 
     72        crs.lonAxis=0 
     73        crs.latAxis=1 
     74        crs.timeAxis=3 
     75         
     76         
    5877        #define unknown 1D CRS: 
    5978        crs=CRSystem(srsName='ndg:crs:unknown', axes=['unknown']) 
    6079        self.systems['ndg:crs:unknown']=crs 
    61          
     80                         
    6281        #define unknown 3D CRS: 
    6382        crs=CRSystem(srsName='ndg:crs:unknown2d', axes=['unknown','unknown']) 
     
    81100     
    82101    def getUnitType(self, unit): 
    83         unittype='unknown' 
     102        #unittype='unknown' 
    84103        if string.lower(unit) in ['second', 'seconds', 's', 'mins','minute','minutes','hour','hours','h','hr','hrs','day','days']: 
    85104            unittype='time' 
    86105        elif string.lower(unit)[:10] in ['days since', 'seconds si', 'minutes si', 'hours sinc','months sin', 'years sinc']: 
    87106            unittype='time' 
    88         elif string.lower(unit) in ['mbar', 'pa','level']: 
     107        elif string.lower(unit) in ['mbar', 'pa']: 
    89108            unittype='pressure' 
    90109        elif string.lower(unit) in ['m', 'km']: 
     
    94113        elif string.lower(unit) in ['degrees_east','degrees_west']: 
    95114            unittype='longitude' 
     115        else: 
     116            unittype = unit 
    96117        return unittype 
    97118     
     
    107128        if knownCRSAxes is not None: 
    108129            axes=knownCRSAxes 
    109                      
    110         if len(axes)==2: 
    111             crs=self.systems['ndg:crs:unknown2d'] 
    112         if len(axes)==3: 
    113             crs=self.systems['ndg:crs:unknown3d'] 
    114         elif len(axes)==4: 
    115             crs=self.systems['ndg:crs:unknown4d'] 
    116         else : 
    117             crs=self.systems['ndg:crs:unknown'] 
     130        crs =None 
     131        #if len(axes)==2: 
     132            #crs=self.systems['ndg:crs:unknown2d'] 
     133        #if len(axes)==3: 
     134            #crs=self.systems['ndg:crs:unknown3d'] 
     135        #elif len(axes)==4: 
     136            #crs=self.systems['ndg:crs:unknown4d'] 
     137        #else : 
     138            #crs=self.systems['ndg:crs:unknown'] 
    118139         
    119140        #this can be extended to accomodate more units        
     
    126147                unittype=self.getUnitType(unit)             
    127148                crsMap.append(unittype) 
    128          
     149            print '-------' 
     150            print axes 
     151            print crsMap 
     152            print '-------' 
     153             
    129154         
    130155        #now try and match up crsMap with known crsystems.         
     
    147172                crs=self.systems[system] 
    148173                break 
    149                          
     174        if match==1:  #couldn't find a match - generate a psuedo crs     
     175            for psystem in self.pseudosystems: 
     176                if len(self.pseudosystems[psystem].axes) != len(crsMap): 
     177                    #wrong dimensionality, skip this crs 
     178                    continue 
     179                axisorder={} 
     180                indexes=[] 
     181                for ax in crsMap: 
     182                    if ax in self.pseudosystems[psystem].axes: 
     183                        axisorder[ax]=crsMap.index(ax) 
     184                        indexes.append(crsMap.index(ax)) 
     185                print 'axorder %s'%axisorder 
     186                print 'crsmap %s'%crsMap 
     187                print 'system %s'%self.pseudosystems[psystem].axes 
     188                ##############################this bit is broken############################ 
     189                if len(axisorder)== len(self.pseudosystems[psystem].axes) - 1: 
     190                    for ax in crsMap: 
     191                        if ax not in axisorder: 
     192                            for i in range(3): 
     193                                if i not in indexes: 
     194                                    axisorder[ax]=i 
     195                                    break           
     196                    
     197                    newaxes={} 
     198                    axisLabels='' 
     199                    print 'axisorder now %s'%axisorder 
     200                    for ax in axisorder: 
     201                        if ax in crsMap: 
     202                            print 'yes' 
     203                            unitpos = crsMap.index(ax) 
     204                            axisname=axes[unitpos] 
     205                            newaxes[axisname]=unitpos 
     206                    print 'newaxes %s'%newaxes 
     207                    for ax in newaxes: 
     208                        axisLabels=axisLabels + ax + ' '        
     209                 #######################################################################    
     210                    crs=CRSystem('ndg:psuedoCRS', axes) 
     211                    crs.axisLabels=axisLabels 
     212                    crs.timeAxis=self.pseudosystems[psystem].timeAxis 
     213                    crs.latAxis=self.pseudosystems[psystem].latAxis 
     214                    crs.latAxis=self.pseudosystems[psystem].lonAxis 
     215        print 'labels %s'%crs.axisLabels 
    150216        return crs, axisorder 
    151217         
    152218def main(): 
    153219    import sys 
    154     sys.stdout = open("/home/dom/Desktop/SVN/trunk/csml/testfiles/gridseries/crsout.txt", "w") 
     220    #sys.stdout = open("/home/dom/svn/trunk/csml/testfiles/gridseries/crsout.txt", "w") 
    155221    cat=CRSCatalogue() 
    156222   
    157223    #test every crs in this netcdf file 
    158224    crsdictionary={} 
    159     f=cdms.open('/home/dom/Desktop/SVN/trunk/csml/testfiles/gridseries/xaaqda@pxs19c1.nc') 
     225    f=cdms.open('/home/dom/svn/trunk/csml/testfiles/gridseries/xaaqda@pxs09c1.nc') 
    160226    for var in f.variables: 
    161227        v=f.variables[var] 
     
    166232            units.append(ax.attributes['units']) 
    167233        crsdictionary[str(axes)]=units 
    168     n=2 
     234    n=0 
    169235    #test each definition in dictionary 
    170236    for key in crsdictionary: 
     
    180246        print 'UNITS: %s'%uns 
    181247        crs, axisorder= cat.determineCRS(axes=axs, units=uns) 
    182         print 'CRS: %s'%crs.srsName 
     248        #print 'CRS: %s'%crs.srsName 
    183249        print 'ORDER:  %s'%axisorder 
    184250        #print 'axis order =%s'%axisorder 
     
    199265        #print '<axisLabels>%s</axisLabels>'%labels 
    200266        #for axis in axes: 
    201             #print '<GridOrdinate>' 
     267            #print '<GridOrdinate>'I just had this same problem with a Kubuntu/Xubuntu dual boot and this solved it for  
    202268  
    203269            #unit=units[axes.index(axis)] 
Note: See TracChangeset for help on using the changeset viewer.