Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
''' Create VRT
Parameters
-----------
fileName : string
gdalDataset : gdal dataset
gdalMetadata : gdal metadata
latlonGrid : numpy 2 layered 2D array with lat/lons of desired grid
'''
# test if input files is ASCAT
iDir, iFile = os.path.split(fileName)
iFileName, iFileExt = os.path.splitext(iFile)
try:
assert iFileName[0:6] == 'ascat_' and iFileExt == '.nc'
except:
raise WrongMapperError
# Create geolocation
subDataset = gdal.Open('NETCDF:"' + fileName + '":lat')
self.GeolocVRT = VRT(srcRasterXSize=subDataset.RasterXSize,
srcRasterYSize=subDataset.RasterYSize)
GeolocMetaDict = [{'src': {'SourceFilename': ('NETCDF:"' + fileName +
'":lon'),
'SourceBand': 1,
'ScaleRatio': 0.00001,
'ScaleOffset': -360},
'dst': {}},
{'src': {'SourceFilename': ('NETCDF:"' + fileName +
'":lat'),
'SourceBand': 1,
'ScaleRatio': 0.00001,
def __init__(self, fileName, gdalDataset, gdalMetadata, **kwargs):
''' Create CSKS VRT '''
if fileName.split('/')[-1][0:4] != "CSKS":
raise WrongMapperError
# Get coordinates
metadata = gdalMetadata['Estimated_Bottom_Left_Geodetic_Coordinates']
bottom_left_lon = float(metadata.split(' ')[1])
bottom_left_lat = float(metadata.split(' ')[0])
metadata = gdalMetadata['Estimated_Bottom_Right_Geodetic_Coordinates']
bottom_right_lon = float(metadata.split(' ')[1])
bottom_right_lat = float(metadata.split(' ')[0])
metadata = gdalMetadata['Estimated_Top_Left_Geodetic_Coordinates']
top_left_lon = float(metadata.split(' ')[1])
top_left_lat = float(metadata.split(' ')[0])
metadata = gdalMetadata['Estimated_Top_Right_Geodetic_Coordinates']
top_right_lon = float(metadata.split(' ')[1])
top_right_lat = float(metadata.split(' ')[0])
metadata = gdalMetadata['Scene_Centre_Geodetic_Coordinates']
center_lon = float(metadata.split(' ')[1])
if (not gdalMetadata or
not 'SATELLITE_IDENTIFIER' in gdalMetadata.keys()):
raise WrongMapperError(filename)
elif gdalMetadata['SATELLITE_IDENTIFIER'] != 'RADARSAT-2':
raise WrongMapperError(filename)
if zipfile.is_zipfile(inputFileName):
# Open product.xml to get additional metadata
zz = zipfile.ZipFile(inputFileName)
productXmlName = os.path.join(os.path.basename(inputFileName).split('.')[0],'product.xml')
productXml = zz.open(productXmlName).read()
else:
# product.xml to get additionali metadata
productXmlName = os.path.join(filename,'product.xml')
if not os.path.isfile(productXmlName):
raise WrongMapperError(filename)
productXml = open(productXmlName).read()
if not IMPORT_SCIPY:
raise NansatReadError(' Radarsat-2 data cannot be read because scipy is not installed! '
' Please do: conda -c conda-forge install scipy ')
# parse product.XML
rs2_0 = Node.create(productXml)
if xmlonly:
self.init_from_xml(rs2_0)
return
# Get additional metadata from product.xml
rs2_1 = rs2_0.node('sourceAttributes')
rs2_2 = rs2_1.node('radarParameters')
def __init__(self, fileName, gdalDataset, gdalMetadata, **kwargs):
''' Create MODIS_L1 VRT '''
# raise error in case of not GOCI L1B
try:
title = gdalMetadata['HDFEOS_POINTS_Scene_Header_Scene_Title']
except (TypeError, KeyError):
raise WrongMapperError
if not title == 'GOCI Level-1B Data':
raise WrongMapperError
# set GOCI projection parameters
lat_0 = gdalMetadata['HDFEOS_POINTS_Map_Projection_Central_Latitude_(parallel)']
lon_0 = gdalMetadata['HDFEOS_POINTS_Map_Projection_Central_Longitude_(meridian)']
rasterXSize = int(gdalMetadata['HDFEOS_POINTS_Scene_Header_number_of_columns'].split(' ')[0])
rasterYSize = int(gdalMetadata['HDFEOS_POINTS_Scene_Header_number_of_rows'].split(' ')[0])
proj4 = '+proj=ortho +lat_0=%s +lon_0=%s units=m +ellps=WGS84 +datum=WGS84 +no_defs' % (lat_0, lon_0)
srs = osr.SpatialReference()
srs.ImportFromProj4(proj4)
projection = srs.ExportToWkt()
geoTransform = (-1391500.0, 500.0, 0.0, 1349500.0, 0.0, -500.0)
# create empty VRT dataset with georeference only
VRT.__init__(self, srcGeoTransform=geoTransform,
srcProjection=projection,
srcRasterXSize=rasterXSize,
spectra (OSW), ocean wind field (OWI), or radial surface velocity
(RVL) (RVL is the default)
GCP_COUNT : int
number of GCPs along each dimention
'''
fPathName, fExt = os.path.splitext(fileName)
# List of Sentinel-1 level-2 components
unwanted_product_components = ['osw', 'owi', 'rvl']
# Remove requested 'product_type' from list of unwanted
unwanted_product_components.pop(unwanted_product_components.index(
product_type.lower()))
# Check if it is Sentinel-1 (or ASAR) level-2 (in S1 data format)
if not gdalMetadata or not 'NC_GLOBAL' in gdalMetadata.keys():
raise WrongMapperError
else:
title = gdalMetadata['NC_GLOBAL#TITLE']
# Raise error if it is not Sentinel-1 format
if not 'Sentinel-1' or 'ASA' in title:
raise WrongMapperError
metadata = {}
for key, val in gdalMetadata.iteritems():
new_key = key.split('#')[-1]
metadata[new_key] = val
subDatasets = gdalDataset.GetSubDatasets()
fileNames = [f[0] for f in subDatasets]
rm_bands = []
def __init__(self, inputFileName, gdalDataset, gdalMetadata,
xmlonly=False, **kwargs):
''' Create Radarsat2 VRT '''
fPathName, fExt = os.path.splitext(inputFileName)
if zipfile.is_zipfile(inputFileName):
# Open zip file using VSI
fPath, fName = os.path.split(fPathName)
filename = '/vsizip/%s/%s' % (inputFileName, fName)
if not 'RS' in fName[0:2]:
raise WrongMapperError('%s: Provided data is not Radarsat-2'
%fName)
gdalDataset = gdal.Open(filename)
gdalMetadata = gdalDataset.GetMetadata()
else:
filename = inputFileName
#if it is not RADARSAT-2, return
if (not gdalMetadata or
not 'SATELLITE_IDENTIFIER' in gdalMetadata.keys()):
raise WrongMapperError(filename)
elif gdalMetadata['SATELLITE_IDENTIFIER'] != 'RADARSAT-2':
raise WrongMapperError(filename)
if zipfile.is_zipfile(inputFileName):
# Open product.xml to get additional metadata
zz = zipfile.ZipFile(inputFileName)
gdalDatasetTmp = gdal.Open(sourceFilename)
# keep name, GDALDataset and size
bandFileNames.append(sourceFilename)
bandSizes.append(gdalDatasetTmp.RasterXSize)
bandDatasets.append(gdalDatasetTmp)
# get mtl file
mtlFiles = glob.glob(coreName+'*[mM][tT][lL].[tT][xX][tT]')
if len(mtlFiles) > 0:
mtlFileName = mtlFiles[0]
else:
raise WrongMapperError
# if not TIF files found - not appropriate mapper
if not bandFileNames:
raise WrongMapperError
# get appropriate band size based on number of unique size and
# required resoltuion
if resolution == 'low':
bandXSise = min(bandSizes)
elif resolution in ['high', 'hi']:
bandXSise = max(bandSizes)
else:
raise OptionError('Wrong resolution %s for file %s' % (resolution, fileName))
# find bands with appropriate size and put to metaDict
metaDict = []
for bandFileName, bandSize, bandDataset in zip(bandFileNames,
bandSizes,
bandDatasets):
if bandSize == bandXSise:
if iMapper == 'mapper_generic' and len(importErrors) > 0:
self.logger.error('\nWarning! '
'The following mappers failed:')
for ie in importErrors:
self.logger.error(importErrors)
# create a Mapper object and get VRT dataset from it
try:
tmpVRT = nansatMappers[iMapper](self.fileName,
gdalDataset,
metadata,
**kwargs)
self.logger.info('Mapper %s - success!' % iMapper)
self.mapper = iMapper.replace('mapper_', '')
break
except WrongMapperError:
pass
# if no mapper fits, make simple copy of the input DS into a VSI/VRT
if tmpVRT is None and gdalDataset is not None:
self.logger.warning('No mapper fits, returning GDAL bands!')
tmpVRT = VRT(gdalDataset=gdalDataset)
for iBand in range(gdalDataset.RasterCount):
tmpVRT._create_band({'SourceFilename': self.fileName,
'SourceBand': iBand + 1})
tmpVRT.dataset.FlushCache()
# if GDAL cannot open the file, and no mappers exist which can make VRT
if tmpVRT is None and gdalDataset is None:
# check if given data file exists
if not os.path.isfile(self.fileName):
raise IOError('%s: File does not exist' % (self.fileName))
def __init__(self, fileName, gdalDataset, gdalMetadata, **kwargs):
''' Create VRT '''
fileBaseName = os.path.basename(fileName)
if not fileBaseName == 'MOD44W.vrt':
raise WrongMapperError
metaDict = [{'src': {'SourceFilename': fileName, 'SourceBand': 1},
'dst': {'wkv': 'land_binary_mask'}}]
# create empty VRT dataset with geolocation only
VRT.__init__(self, gdalDataset)
# add bands with metadata and corresponding values to the empty VRT
self._create_bands(metaDict)
mm = pti.get_gcmd_instrument('MODIS')
ee = pti.get_gcmd_platform('TERRA')
self.dataset.SetMetadataItem('instrument', json.dumps(mm))
self.dataset.SetMetadataItem('platform', json.dumps(ee))
def __init__(self, fileName, gdalDataset, gdalMetadata, logLevel=30,
**kwargs):
if fileName[0:len(keywordBase)] != keywordBase:
raise WrongMapperError(__file__,
"Not Nora10 data converted from felt to netCDF")
requestedTime = datetime.strptime(fileName[len(keywordBase)+1:],
'%Y%m%d%H%M')
# For correct rounding
fileTime = requestedTime + timedelta(minutes=30)
fileTime = fileTime - timedelta(minutes=fileTime.minute)
nc_file = (baseFolder + 'windspeed_10m' +
fileTime.strftime('/%Y/%m/') + 'windspeed_' +
fileTime.strftime('%Y%m%d%H.nc'))
nc_file_winddir = (baseFolder + 'winddir_10m' +
fileTime.strftime('/%Y/%m/') +
'winddir_' + fileTime.strftime('%Y%m%d%H.nc'))
# Would prefer to use geotransform, but ob_tran