Python Nansat.watermask方法代码示例

# 需要导入模块: from nansat import Nansat [as 别名]
# 或者: from nansat.Nansat import watermask [as 别名]
 def test_watermask(self):
     ''' if watermask data exists: should fetch array with watermask
         else:                     should raise an error'''
     n1 = Nansat(self.test_file_gcps, logLevel=40)
     mod44path = os.getenv('MOD44WPATH')
     if mod44path is not None and os.path.exists(mod44path + '/MOD44W.vrt'):
         wm = n1.watermask()[1]
         self.assertEqual(type(wm), np.ndarray)
         self.assertEqual(wm.shape[0], n1.shape()[0])
         self.assertEqual(wm.shape[1], n1.shape()[1])

# 需要导入模块: from nansat import Nansat [as 别名]
# 或者: from nansat.Nansat import watermask [as 别名]
# 5. Write image
n.reproject() # 1.
lons, lats = n.get_corners() # 2.
pxlRes = distancelib.getPixelResolution(array(lats), array(lons), n[1])
pxlRes = array(pxlRes)*360/40000 # great circle distance
srsString = "+proj=latlong +datum=WGS84 +ellps=WGS84 +no_defs"
#~ extentString = '-lle %f %f %f %f -ts 3000 3000' % (min(lons), min(lats), max(lons), max(lats))
extentString = '-lle %f %f %f %f -tr %f %f' % (min(lons), min(lats), \
                max(lons), max(lats), pxlRes[1], pxlRes[0])
d = Domain(srs=srsString, ext=extentString) # 3.
print d
n.reproject(d) # 4.

# get array with watermask (landmask) b 
# it must be done after reprojection!
# 1. Get Nansat object with watermask
# 2. Get array from Nansat object. 0 - land, 1 - water
#wm = n.watermask(mod44path='/media/magDesk/media/SOLabNFS/store/auxdata/coastline/mod44w/')
wm = n.watermask(mod44path='/media/data/data/auxdata/coastline/mod44w/')
wmArray = wm[1]

figureName = oPath + fileName + '_proj.png'
n.write_figure(fileName=figureName, clim=[3,133], \
               mask_array=wmArray, mask_lut={0: [204, 153, 25]}) # 5.

#~ # make KML file with image borders (to be opened in Googe Earth)
#~ n.write_kml(kmlFileName=oPath + fileName + '_preview.kml')

# make KML image file with image borders (to be opened in Googe Earth)
n.write_image_kml(kmlFileName=oPath + fileName + '.kml', kmlFigureName=figureName)

# 需要导入模块: from nansat import Nansat [as 别名]
# 或者: from nansat.Nansat import watermask [as 别名]
fig.save(oFileName + '04_logo.png')

# Create a Figure object (fig)
fig = Figure(array)
# Get lat/lon arrays from Nansat object (may take some time)
lonGrid, latGrid = n.get_geolocation_grids()
# Make figure with lat/lon grids
fig.process(cmin=10, cmax=60, latGrid=latGrid, lonGrid=lonGrid,
            latlonGridSpacing=10, latlonLabels=10)
# save the fig
fig.save(oFileName + '05_latlon.png', )

# Create a Figure object (fig)
fig = Figure(array)
# Get Nansat object with watermask
wm = n.watermask()
# Get array from Nansat object. 0 - land, 1 - water
wmArray = wm[1]
# Compute min and max valuse from ratio
clim = fig.clim_from_histogram(ratio=1.0)
# Make figure with land overlay (gray) and apply brightness gamma correction
fig.process(cmin=clim[0], cmax=clim[1], mask_array=wmArray,
               mask_lut={2: [128, 128, 128]}, logarithm=True, gamma=3)
# save the fig
fig.save(oFileName + '06_land.png', )

# create 3D numpy array
array = np.array([n[1], n[2], n[3]])
# Create a Figure object (fig) from 3D array
fig = Figure(array)
# Compute min and max valuse from ratio

# 需要导入模块: from nansat import Nansat [as 别名]
# 或者: from nansat.Nansat import watermask [as 别名]

#.........这里部分代码省略.........
            # try to create Nansat object
            try:
                n = Nansat(iPath + dirName + '/' + fileName, mapperName='asar', logLevel=27)
            except Exception as e:
                print "Failed to create Nansat object:"
                print str(e)
                os.rmdir(oPath + fileName[0:27] + '/' )
                continue
                

            #~ Get the bands
            raw_counts = n[1]
            inc_angle = n[2]

            #~ NICE image (roughness)
            pol = n.bands()[3]['polarization']
            if pol == 'HH':
                ph = (2.20495, -14.3561e-2, 11.28e-4)
                sigma0_hh_ref = exp( ( ph[0]+inc_angle*ph[1]+inc_angle**2*ph[2])*log(10) )
                roughness = n[3]/sigma0_hh_ref
            elif pol == 'VV':
                pv = (2.29373, -15.393e-2, 15.1762e-4)
                sigma0_vv_ref = exp( ( pv[0]+inc_angle*pv[1]+inc_angle**2*pv[2])*log(10) )
                roughness = n[3]/sigma0_vv_ref

            #~ Create new band
            n.add_band(bandID=4, array=roughness, \
               parameters={'name':'roughness', \
               'wkv': 'surface_backwards_scattering_coefficient_of_radar_wave', \
               'dataType': 6})

            # Reproject image into Lat/Lon WGS84 (Simple Cylindrical) projection
            # 1. Cancel previous reprojection
            # 2. Get corners of the image and the pixel resolution
            # 3. Create Domain with stereographic projection, corner coordinates 1000m
            # 4. Reproject
            # 5. Write image
            n.reproject() # 1.
            lons, lats = n.get_corners() # 2.
            # Pixel resolution
            #~ pxlRes = distancelib.getPixelResolution(array(lats), array(lons), n.shape())
            #~ pxlRes = array(pxlRes)*360/40000 # great circle distance
            pxlRes = array(distancelib.getPixelResolution(array(lats), array(lons), n.shape(), 'deg'))
            
            
            ipdb.set_trace()
            
            
            if min(lats) >= 65 and max(lats) >= 75 and max(lats)-min(lats) >= 13:
               pxlRes = array([0.00065, 0.00065])*2 # make the resolution 150x150m
            #~ pxlRes = pxlRes*7 # make the resolution worser
            srsString = "+proj=latlong +datum=WGS84 +ellps=WGS84 +no_defs"
            #~ extentString = '-lle %f %f %f %f -ts 3000 3000' % (min(lons), min(lats), max(lons), max(lats))
            extentString = '-lle %f %f %f %f -tr %f %f' % (min(lons), min(lats), \
                            max(lons), max(lats), pxlRes[1], pxlRes[0])
            d = Domain(srs=srsString, ext=extentString) # 3.
            n.reproject(d) # 4.

            if useMask:
               # get array with watermask (landmask) b
               # it must be done after reprojection!
               # 1. Get Nansat object with watermask
               # 2. Get array from Nansat object. 0 - land, 1 - water
               #wm = n.watermask(mod44path='/media/magDesk/media/SOLabNFS/store/auxdata/coastline/mod44w/')
               wm = n.watermask(mod44path='/media/data/data/auxdata/coastline/mod44w/')
               wmArray = wm[1]

               #~ ОШИБКА numOfColor=255 не маскирует, потому что в figure.apply_mask: availIndeces = range(self.d['numOfColor'], 255 - 1)
               #~ n.write_figure(fileName=figureName, bands=[3], \
                            #~ numOfColor=255, mask_array=wmArray, mask_lut={0: 0},
                            #~ clim=[0,0.15], cmapName='gray', transparency=0) # 5.
               n.write_figure(fileName=figureName, bands=[4], \
                                 mask_array=wmArray, mask_lut={0: [0,0,0]},
                                 clim=[0,2], cmapName='gray', transparency=[0,0,0]) # 5.
            else:
               n.write_figure(fileName=figureName, bands=[1], \
                              clim=[0,2], cmapName='gray', transparency=[0,0,0]) # 5.

            # open the input image and convert to RGBA for further tiling with slbtiles
            input_img = Image.open(figureName)
            output_img = input_img.convert("RGBA")
            output_img.save(figureName)

            # make KML image
            n.write_kml_image(kmlFileName=kmlName, kmlFigureName=figureName)

            #~ Change the file permissions
            os.chmod(oPath, 0777)
            os.chmod(oPath + fileName[0:27] + '/', 0777)
            os.chmod(kmlName, 0777)
            os.chmod(figureName, 0777)

            #~ Change the owner and group
            #~ os.chown(oPath, 1111, 1111)
            #~ os.chown(oPath + fileName[0:27] + '/', 1111, 1111)
            #~ os.chown(kmlName, 1111, 1111)
            #~ os.chown(figureName, 1111, 1111)
            
            #~ garbage collection
            gc.collect()