generic

basic reading of Sentinel-2 data

s2d2.read_sentinel2.dn2toa_s2(I)[source]

convert the digital numbers of Sentinel-2 to top of atmosphere (TOA), see for more details [wwwS2L1C].

Parameters:: I (numpy.ndarray, dim={2,3}, size=(m,n)) – grid with intensities
Returns:: grid with top of atmosphere reflectances
Return type:: numpy.ndarray, dim={2,3}, size=(m,n)

Notes

[wwwS2L1C]

https://sentinel.esa.int/web/sentinel/technical-guides/sentinel-2-msi/level-1c/algorithm

s2d2.read_sentinel2.get_flight_bearing_from_detector_mask_s2(Det)[source]

Parameters:: Det (numpy.array, size=(m,n), ndim={2,3}, dtype=integer) – array with numbers of the different detectors in Sentinel-2
Returns:: ψ – general bearing of the satellite
Return type:: float, unit=degrees

s2d2.read_sentinel2.get_integration_and_sampling_time_s2(ds_path, fname='MTD_DS.xml', s2_dict=None)[source]

Parameters:

ds_path (string) – location where metadata of datastrip is situated
fname (string, default='MTD_DS.xml') – metadata filename
s2_dict (dictionary, default=None) – metadata dictionary, can be used instead of the string based method

Notes

The metadata is scattered over the file structure of Sentinel-2, L1C

* S2X_MSIL1C_20XX...
├ AUX_DATA
│  └ DS_XXX_XXXX...
│     └ QI_DATA
│        └ MTD_DS.xml <- metadata about the data-strip
├ GRANULE
│  └ L1C_TXXXX_XXXX...
│     ├ AUX_DATA
│     ├ IMG_DATA
│     ├ QI_DATA
│     └ MTD_TL.xml <- metadata about the tile
├ HTML
├ rep_info
├ manifest.safe
├ INSPIRE.xml
└ MTD_MSIL1C.xml <- metadata about the product

The following acronyms are used:

DS : datastrip
TL : tile
QI : quality information
AUX : auxiliary
MTD : metadata
MSI : multi spectral instrument
L1C : product specification,i.e.: level 1, processing step C

s2d2.read_sentinel2.list_central_wavelength_msi()[source]

create dataframe with metadata about Sentinel-2

Returns:

metadata and general multispectral information about the MSI instrument that is onboard Sentinel-2, having the following collumns:

center_wavelength, unit=µm : central wavelength of the band

full_width_half_max, unit=µm : extent of the spectral sensativity

bandid : number for identification in the meta data

resolution, unit=m : spatial resolution of a pixel

along_pixel_size, unit=µm : physical size of the sensor

across_pixel_size, unit=µm : size of the photosensative sensor

focal_length, unit=m : lens characteristic of the telescope

field_of_view, unit=degrees : angle of swath of instrument

crossdetector_parallax, unit=degress : in along-track direction

name : general name of the band, if applicable

solar_illumination, unit=W m-2 µm-1 :

mostly following the naming convension of the STAC EO extension [stac]

Return type:

pandas.dataframe

Notes

The Multi Spectral instrument (MSI) has a Tri Mirror Anastigmat (TMA) telescope, which is opened at F/4 and has a focal length of 60 centimeter.

The crossdetector parallex is given here as well, see [La15]. While the crossband parallax is 0.018 for VNIR and 0.010 for the SWIR detector, respectively. The dimensions of the Silicon CMOS detector are given in [MG10]. While the SWIR detector is based on a MCT sensor.

The following acronyms are used:

MSI : multi-spectral instrument

MCT : mercury cadmium telluride, HgCdTe

TMA : trim mirror anastigmat

VNIR : very near infra-red

SWIR : short wave infra-red

CMOS : silicon complementary metal oxide semiconductor, Si

Examples

>>> from s2d2.read_sentinel2 import list_central_wavelength_msi

make a selection by name:

>>> boi = ['red', 'green', 'blue', 'near infrared']
>>> s2_df = list_central_wavelength_msi()
>>> s2_df = s2_df[s2_df['common_name'].isin(boi)]
>>> s2_df
         wavelength  bandwidth  resolution       name  bandid
B02         492         66          10           blue       1
B03         560         36          10          green       2
B04         665         31          10            red       3
B08         833        106          10  near infrared       7

similarly you can also select by pixel resolution:

>>> s2_df = list_central_wavelength_msi()
>>> tw_df = s2_df[s2_df['resolution']==20]
>>> tw_df.index
Index(['B05', 'B06', 'B07', 'B8A', 'B11', 'B12'], dtype='object')

References

[La15]

Languille et al. “Sentinel-2 geometric image quality commissioning: first results” Proceedings of the SPIE, 2015.

[MG10]

Martin-Gonthier et al. “CMOS detectors for space applications: From R&D to operational program with large volume foundry”, Proceedings of the SPIE conference on sensors, systems, and next generation satellites XIV, 2010.

[stac]

https://github.com/stac-extensions/eo

s2d2.read_sentinel2.read_band_s2(path, band=None)[source]

This function takes as input the Sentinel-2 band name and the path of the folder that the images are stored, reads the image and returns the data as an array

Parameters:

path (string) – path of the folder, or full path with filename as well
band (string, optional) – Sentinel-2 band name, for example ‘04’, ‘8A’.

Returns:

data (numpy.array, size=(_,_)) – array of the band image
spatialRef (string) – projection
geoTransform (tuple) – affine transformation coefficients
targetprj (osr.SpatialReference object) – spatial reference

See also

list_central_wavelength_msi: creates a dataframe for the MSI instrument
read_stack_s2: reading several Sentinel-2 bands at once into a stack
s2d2.mapping_input.read_geo_image: basic function to import geographic imagery data

Examples

>>> from s2d2.read_sentinel2 import read_band_s2

>>> path = '/GRANULE/L1C_T15MXV_A027450_20200923T163313/IMG_DATA/'
>>> band = '02'
>>> _,spatialRef,geoTransform,targetprj = read_band_s2(path, band)

>>> spatialRef
'PROJCS["WGS 84 / UTM zone 15S",GEOGCS["WGS 84",DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],
UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],
AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-93],
PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],
PARAMETER["false_northing",10000000],UNIT["metre",1,
AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],
AUTHORITY["EPSG","32715"]]'
>>> geoTransform
(600000.0, 10.0, 0.0, 10000000.0, 0.0, -10.0)
>>> targetprj
<osgeo.osr.SpatialReference; proxy of <Swig Object of type
'OSRSpatialReferenceShadow *' at 0x7f9a63ffe450> >

s2d2.read_sentinel2.read_cloud_mask(path_meta, geoTransform)[source]

Parameters:

path_meta (string) – directory where meta data is situated, ‘MSK_CLOUDS_B00.gml’ is typically the file of interest
geoTransform (tuple) – affine transformation coefficients

Returns:

array which highlights where clouds could be

Return type:

numpy.ndarray, size=(m,n), dtype=int

s2d2.read_sentinel2.read_detector_mask(path_meta, boi, geoTransform)[source]

create array of with detector identification

Sentinel-2 records in a pushbroom fasion, collecting reflectance with a ground resolution of more than 270 kilometers. This data is stacked in the flight direction, and then cut into granules. However, the sensorstrip inside the Multi Spectral Imager (MSI) is composed of several CCD arrays.

This function collects the geometry of these sensor arrays from the meta- data. Since this is stored in a gml-file.

Parameters:

path_meta (string) – path where the meta-data is situated.
boi (pandas.DataFrame) – list with bands of interest
geoTransform (tuple, size={(6,),(8,)}) – affine transformation coefficients

Returns:

numpy.ma.array, size=(msk_dim[0 – array where each pixel has the ID of the detector, of a specific band

Return type:

1],len(boi)), dtype=int8

See also

list_central_wavelength_msi: creates a dataframe for the MSI instrument
read_view_angles_s2: read grid of Sentinel-2 observation angles

Notes

The metadata is scattered over the file structure of Sentinel-2, L1C

* S2X_MSIL1C_20XX...
├ AUX_DATA
├ DATASTRIP
│  └ DS_XXX_XXXX...
│     └ QI_DATA
│        └ MTD_DS.xml <- metadata about the data-strip
├ GRANULE
│  └ L1C_TXXXX_XXXX...
│     ├ AUX_DATA
│     ├ IMG_DATA
│     ├ QI_DATA
│     └ MTD_TL.xml <- metadata about the tile
├ HTML
├ rep_info
├ manifest.safe
├ INSPIRE.xml
└ MTD_MSIL1C.xml <- metadata about the product

The following acronyms are used:

DS : datastrip
TL : tile
QI : quality information
AUX : auxiliary
MTD : metadata
MSI : multi spectral instrument
L1C : product specification,i.e.: level 1, processing step C

Examples

>>> from s2d2.read_sentinel2 import (
>>>    list_central_wavelength_msi, read_detector_mask
>>> )

>>> path_meta = '/GRANULE/L1C_T15MXV_A027450_20200923T163313/QI_DATA'
>>> boi = ['red', 'green', 'blue', 'near infrared']
>>> s2_df = list_central_wavelength_msi()
>>> boi_df = s2_df[s2_df['common_name'].isin(boi)]
>>> geoTransform = (600000.0, 10.0, 0.0, 10000000.0, 0.0, -10.0)
>>>
>>> det_stack = read_detector_mask(path_meta, boi_df, geoTransform)

s2d2.read_sentinel2.read_detector_time_s2(path, fname='MTD_DS.xml', s2_df=None)[source]

get the detector metadata of the relative detector timing

Parameters:

path (string) – path where the meta-data is situated
fname (string) – file name of the metadata.

Returns:

det_time (numpy.array, numpy.datetime64) – detector time
det_name (list, size=(13,)) – list of the different detector codes
det_meta (numpy.array, size=(13,4)) –
metadata of the individual detectors, each detector has the following:
1. spatial resolution [m]
2. minimal spectral range [nm]
3. maximum spectral range [nm]
4. mean spectral range [nm]
line_period (float, numpy.timedelta64) – temporal sampling distance

Notes

The sensor blocks are arranged as follows, with ~98 pixels overlap:

    ┌-----┐   ┌-----┐   ┌-----┐   ┌-----┐   ┌-----┐   ┌-----┐
    |DET02|   |DET04|   |DET06|   |SCA08|   |SCA10|   |SCA12|
    └-----┘   └-----┘   └-----┘   └-----┘   └-----┘   └-----┘
┌-----┐   ┌-----┐   ┌-----┐   ┌-----┐   ┌-----┐   ┌-----┐
|DET01|   |DET03|   |SCA05|   |SCA07|   |SCA09|   |SCA11|
└-----┘   └-----┘   └-----┘   └-----┘   └-----┘   └-----┘

a forward looking array, while band order is reverse for aft looking
<-2592 10m pixels->
┌-----------------┐    #*# satellite
|       B02       |     |
├-----------------┤     | flight
|       B08       |     | direction
├-----------------┤     |
|       B03       |     v
├-----------------┤
       etc.
the detector order is B02, B08, B03, B10, B04, B05,
B11, B06, B07, B8A, B12, B01 and B09

s2d2.read_sentinel2.read_geotransform_s2(path, fname='MTD_TL.xml', resolution=10)[source]

get the mapping transformation of the Sentinel-2 image

Parameters:

path (string) – location where the meta data is situated
fname (string) – file name of the meta-data file
resolution ({float,integer}, unit=meters, default=10) – resolution of the grid

Returns:

geoTransform – affine transformation coefficients

Return type:

tuple, size=(1,6)

See also

s2d2.image_coordinate_tools.map2pix: transform map- to image-coordinates
s2d2.image_coordinate_tools.map2pix: transform image- to map-coordinates

Notes

The metadata is scattered over the file structure of Sentinel-2, L1C

* S2X_MSIL1C_20XX...
├ AUX_DATA
├ DATASTRIP
│  └ DS_XXX_XXXX...
│     └ QI_DATA
│        └ MTD_DS.xml <- metadata about the data-strip
├ GRANULE
│  └ L1C_TXXXX_XXXX...
│     ├ AUX_DATA
│     ├ IMG_DATA
│     ├ QI_DATA
│     └ MTD_TL.xml <- metadata about the tile
├ HTML
├ rep_info
├ manifest.safe
├ INSPIRE.xml
└ MTD_MSIL1C.xml <- metadata about the product

The metadata structure is as follows:

MTD_TL.xml

└ n1:Level-1C_Tile_ID: ├ n1:General_Info ├ n1:Geometric_Info │ ├ Tile_Geocoding │ │ ├ HORIZONTAL_CS_NAME │ │ ├ HORIZONTAL_CS_CODE │ │ ├ Size : resolution={“10”,”20”,”60”} │ │ │ ├ NROWS │ │ │ └ NCOLS │ │ └ Geoposition │ │ ├ ULX │ │ ├ ULY │ │ ├ XDIM │ │ └ YDIM │ └ Tile_Angles └ n1:Quality_Indicators_Info

The following acronyms are used:

DS : datastrip
TL : tile
QI : quality information
AUX : auxiliary
MTD : metadata
MSI : multi spectral instrument
L1C : product specification,i.e.: level 1, processing step C

s2d2.read_sentinel2.read_mean_sun_angles_s2(path, fname='MTD_TL.xml')[source]

Read the xml-file of the Sentinel-2 scene and extract the mean sun angles.

Parameters:

path (string) – path where xml-file of Sentinel-2 is situated

Returns:

Zn (float) – Mean solar zentih angle of the scene, in degrees.
Az (float) – Mean solar azimuth angle of the scene, in degrees

Notes

The azimuth angle declared in the following coordinate frame:

     ^ North & y
     |
- <--|--> +
     |
     └----> East & x

The angles related to the sun are as follows:

surface normal              * sun
^                     ^    /
|                     |   /
|-- zenith angle      |  /
| /                   | /|
|/                    |/ | elevation angle
└----                 └------

The metadata is scattered over the file structure of Sentinel-2, L1C

* S2X_MSIL1C_20XX...
├ AUX_DATA
├ DATASTRIP
│  └ DS_XXX_XXXX...
│     └ QI_DATA
│        └ MTD_DS.xml <- metadata about the data-strip
├ GRANULE
│  └ L1C_TXXXX_XXXX...
│     ├ AUX_DATA
│     ├ IMG_DATA
│     ├ QI_DATA
│     └ MTD_TL.xml <- metadata about the tile
├ HTML
├ rep_info
├ manifest.safe
├ INSPIRE.xml
└ MTD_MSIL1C.xml <- metadata about the product

The following acronyms are used:

s2 : Sentinel-2
DS : datastrip
TL : tile
QI : quality information
AUX : auxiliary
MTD : metadata
MSI : multi spectral instrument
L1C : product specification,i.e.: level 1, processing step C

s2d2.read_sentinel2.read_orbit_number_s2(path, fname='MTD_MSIL1C.xml')[source]: read the orbit number of the Sentinel-2 image from the metadata

s2d2.read_sentinel2.read_sensing_time_s2(path, fname='MTD_TL.xml')[source]

Parameters:

path (string) – path where the meta-data is situated
fname (string) – file name of the metadata.

Returns:

time of image acquisition

Return type:

numpy.datetime64

See also

get_s2_granual_id

Notes

The recording time is the average sensing within the tile, which is a combination of forward and backward looking datastrips. Schematically, this might look like:

x recording time of datastrip, i.e.: the start of the strip
× recording time of tile, i.e.: the weighted average
                           _________  datastrip
                 ____x____/        /_________
                /        /        /        /
               /        /        /        /
              /        /        /        /
             /        /        /        /
            /        /        /        /
           /        /        /        /
          /    ┌---/------┐ /        /
         /     |  /       |/        /
        /      | /  ×     /        /
       /       |/        /|       /
      /        /--------/-┘      /
     /        /  tile  /        /
    /        /        /        /
   /        /        /        /
  /        /        /        /
 /        /        /        /
/        _________/        /
_________        /________/

Examples

demonstrate the code when in a Scihub data structure

>>> import os
>>> from s2d2.handler.sentinel2 import get_s2_image_locations
>>> from s2d2.read_sentinel2 import read_sensing_time_s2

>>> fpath = '/Users/Data/'
>>> sname = 'S2A_MSIL1C_20200923T163311_N0209_R140_T15MXV_20200923T200821.SAFE'
>>> fname = 'MTD_MSIL1C.xml'
>>> s2_path = os.path.join(fpath, sname, fname)
>>> s2_df,_ = get_s2_image_locations(fname, s2_df)
>>> s2_df['filepath']
B02 'GRANULE/L1C_T15MXV_A027450_20200923T163313/IMG_DATA/T115MXV...'
B03 'GRANULE/L1C_T15MXV_A027450_20200923T163313/IMG_DATA/T15MXV...'
>>> full_path = '/'.join(s2_df['filepath'][0].split('/')[:-2])
'GRANULE/L1C_T15MXV_A027450_20200923T163313'
>>> rec_time = read_sensing_time_s2(path, fname='MTD_TL.xml')
>>> rec_time

s2d2.read_sentinel2.read_stack_s2(s2_df)[source]

read imagery data of interest into an three dimensional np.array

Parameters:

s2_df (pandas.Dataframe) – metadata and general multispectral information about the MSI instrument that is onboard Sentinel-2

Returns:

im_stack (numpy.ndarray, ndim=3) – array of the band image
spatialRef (string) – projection
geoTransform (tuple, size=(8,1)) – affine transformation coefficients
targetprj (osr.SpatialReference object) – spatial reference

See also

list_central_wavelength_msi: creates a dataframe for the MSI instrument
get_s2_image_locations: provides dataframe with specific file locations
read_band_s2: reading a single Sentinel-2 band

Examples

>>> s2_df = list_central_wavelength_msi()
>>> s2_df = s2_df[s2_df['gsd'] == 10]
>>> s2_df,_ = get_s2_image_locations(IM_PATH, s2_df)

>>> im_stack, spatialRef, geoTransform, targetprj = read_stack_s2(s2_df)

s2d2.read_sentinel2.read_sun_angles_s2(path, fname='MTD_TL.xml')[source]

This function reads the xml-file of the Sentinel-2 scene and extracts an array with sun angles, as these vary along the scene.

Parameters:

path (string) – path where xml-file of Sentinel-2 is situated

Returns:

Zn (numpy.array, size=(m,n), dtype=float) – array of the solar zenith angles, in degrees.
Az (numpy.array, size=(m,n), dtype=float) – array of the solar azimuth angles, in degrees.

Notes

The computation of the solar angles is done in two steps: 1) Computation of the solar angles in J2000; 2) Transformation of the vector to the mapping frame.

1. The outputs of the first step is the solar direction normalized vector with the Earth-Sun distance, considering that the direction of the sun is the same at the centre of the Earth and at the centre of the Sentinel-2 satellite.
1. Attitude of the satellite platform are used to rotate the solar vector to the mapping frame. Also Ground Image Calibration Parameters (GICP Diffuser Model) are used to transform from the satellite to the diffuser, as Sentinel-2 has a forward and backward looking sensor configuration. The diffuser model also has stray-light correction and a Bi-Directional Reflection Function model.

The angle(s) are declared in the following coordinate frame:

     ^ North & y
     |
- <--|--> +
     |
     +----> East & x

The angles related to the sun are as follows:

surface normal              * sun
^                     ^    /
|                     |   /
|-- zenith angle      |  /
| /                   | /|
|/                    |/ | elevation angle
+----                 +------

Two different coordinate system are used here:

indexing   |           indexing    ^ y
system 'ij'|           system 'xy' |
           |                       |
           |       i               |       x
   --------+-------->      --------+-------->
           |                       |
           |                       |
image      | j         map         |
based      v           based       |

The metadata is scattered over the file structure of Sentinel-2, L1C

* S2X_MSIL1C_20XX...
├ AUX_DATA
├ DATASTRIP
│  └ DS_XXX_XXXX...
│     └ QI_DATA
│        └ MTD_DS.xml <- metadata about the data-strip
├ GRANULE
│  └ L1C_TXXXX_XXXX...
│     ├ AUX_DATA
│     ├ IMG_DATA
│     ├ QI_DATA
│     └ MTD_TL.xml <- metadata about the tile
├ HTML
├ rep_info
├ manifest.safe
├ INSPIRE.xml
└ MTD_MSIL1C.xml <- metadata about the product

The metadata structure looks like:

* MTD_TL.xml
└ n1:Level-1C_Tile_ID
   ├ n1:General_Info
   ├ n1:Geometric_Info
   │  ├ Tile_Geocoding
   │  └ Tile_Angles
   │     ├ Sun_Angles_Grid
   │     │  ├ Zenith
   │     │  │  ├ COL_STEP
   │     │  │  ├ ROW_STEP
   │     │  │  └ Values_List
   │     │  └ Azimuth
   │     │     ├ COL_STEP
   │     │     ├ ROW_STEP
   │     │     └ Values_List
   │     ├ Mean_Sun_Angle
   │     └ Viewing_Incidence_Angles_Grids
   └ n1:Quality_Indicators_Info

The following acronyms are used:

s2 : Sentinel-2
DS : datastrip
TL : tile
QI : quality information
AUX : auxiliary
MTD : metadata
MSI : multi spectral instrument
L1C : product specification,i.e.: level 1, processing step C

s2d2.read_sentinel2.read_view_angles_s2(path, fname='MTD_TL.xml', det_stack=array([], dtype=float64), boi_df=None)[source]

This function reads the xml-file of the Sentinel-2 scene and extracts an array with viewing angles of the MSI instrument.

Parameters:

path (string) – path where xml-file of Sentinel-2 is situated
fname (string) – the name of the metadata file, sometimes this is changed
boi_df (pandas.dataframe, default=4) – each band has a somewhat minute but different view angle

Returns:

Zn (numpy.ma.array, size=(m,n), dtype=float) – masked array of the solar zenith angles, in degrees.
Az (numpy.ma.array, size=(m,n), dtype=float) – masked array of the solar azimuth angles, in degrees.

See also

list_central_wavelength_s2

Notes

The azimuth angle is declared in the following coordinate frame:

     ^ North & y
     |
- <--|--> +
     |
     └----> East & x

The angles related to the satellite are as follows:

     #*#                   #*# satellite
^    /                ^    /|
|   /                 |   / | nadir
|-- zenith angle      |  /  v
| /                   | /|
|/                    |/ | elevation angle
└----- surface        └------

The metadata is scattered over the file structure of Sentinel-2, L1C

* S2X_MSIL1C_20XX...
├ AUX_DATA
├ DATASTRIP
│  └ DS_XXX_XXXX...
│     └ QI_DATA
│        └ MTD_DS.xml <- metadata about the data-strip
├ GRANULE
│  └ L1C_TXXXX_XXXX...
│     ├ AUX_DATA
│     ├ IMG_DATA
│     ├ QI_DATA
│     └ MTD_TL.xml <- metadata about the tile
├ HTML
├ rep_info
├ manifest.safe
├ INSPIRE.xml
└ MTD_MSIL1C.xml <- metadata about the product

The metadata structure looks like:

* MTD_TL.xml
└ n1:Level-1C_Tile_ID
   ├ n1:General_Info
   ├ n1:Geometric_Info
   │  ├ Tile_Geocoding
   │  └ Tile_Angles
   │     ├ Sun_Angles_Grid
   │     ├ Mean_Sun_Angle
   │     └ Viewing_Incidence_Angles_Grids
   │        ├ Zenith
   │        │  ├ COL_STEP
   │        │  ├ ROW_STEP
   │        │  └ Values_List
   │        └ Azimuth
   │           ├ COL_STEP
   │           ├ ROW_STEP
   │           └ Values_List
   └ n1:Quality_Indicators_Info

The following acronyms are used:

s2 : Sentinel-2
DS : datastrip
TL : tile
QI : quality information
AUX : auxiliary
MTD : metadata
MSI : multi spectral instrument
L1C : product specification,i.e.: level 1, processing step C

advanced reading of Sentinel-2 sensors

s2d2.sensor_readings_sentinel2.get_flight_bearing_from_gnss_s2(path, spatialRef, rec_tim, fname='MTD_DS.xml')[source]

get the direction/argument/heading of the Sentinel-2 acquisition

Parameters:

path (string) – directory where the meta-data is located
spatialRef (osgeo.osr.SpatialReference) – projection system used
rec_tim ({integer,float}) – time stamp of interest
fname (string) – filename of the meta-data file

Returns:

az – array with argument values, based upon the map projection given

Return type:

numpy.array, size=(m,1), float

See also

get_s2_image_locations: to get the datastrip id

Notes

The metadata is scattered over the file structure of Sentinel-2, L1C

* S2X_MSIL1C_20XX...
├ AUX_DATA
├ DATASTRIP
│  └ DS_XXX_XXXX...
│     └ QI_DATA
│        └ MTD_DS.xml <- metadata about the data-strip
├ GRANULE
│  └ L1C_TXXXX_XXXX...
│     ├ AUX_DATA
│     ├ IMG_DATA
│     ├ QI_DATA
│     └ MTD_TL.xml <- metadata about the tile
├ HTML
├ rep_info
├ manifest.safe
├ INSPIRE.xml
└ MTD_MSIL1C.xml <- metadata about the product

The following acronyms are used:

DS : datastrip
TL : tile
QI : quality information
AUX : auxiliary
MTD : metadata
MSI : multi spectral instrument
L1C : product specification,i.e.: level 1, processing step C

s2d2.sensor_readings_sentinel2.get_flight_orientation_s2(ds_path, fname='MTD_DS.xml', s2_dict=None)[source]

get the flight path and orientations of the Sentinel-2 satellite during acquisition.

It is also possible to only give the dictionary, as this has such metadata within.

Parameters:

ds_path (string) – directory where the meta-data is located
fname (string, default='MTD_DS.xml') – filename of the meta-data file
s2_dict (dictonary, default=None) – metadata of the Sentinel-2 platform

Returns:

sat_time (numpy.array, size=(m,1), unit=nanosec) – satellite timestamps
sat_angles (numpy.array, size=(m,3)) – satellite orientation angles, given in ECEF
s2_dict (dictonary) – updated with keys: “time”, “quat”

See also

get_flight_path_s2: get positions of the flight path

s2d2.sensor_readings_sentinel2.get_flight_path_s2(ds_path, fname='MTD_DS.xml', s2_dict=None)[source]

It is also possible to only give the dictionary, as this has such metadata within.

Parameters:

ds_path (string) – location of the metadata file
fname (string, default='MTD_DS.xml') – name of the xml-file that has the metadata
s2_dict (dictonary) – metadata of the Sentinel-2 platform

Returns:

sat_time (numpy.array, size=(m,1), dtype=np.datetime64, unit=ns) – time stamp of the satellite positions
sat_xyz (numpy.array, size=(m,3), dtype=float, unit=meter) – 3D coordinates of the satellite within an Earth centered Earth fixed (ECEF) frame.
sat_err (numpy.array, size=(m,3), dtype=float, unit=meter) – error estimate of the 3D coordinates given by “sat_xyz”
sat_uvw (numpy.array, size=(m,3), dtype=float, unit=meter sec-1) – 3D velocity vectors of the satellite within an Earth centered Earth fixed (ECEF) frame.
s2_dict (dictonary) – updated with keys: “gps_xyz”, “gps_uvw”, “gps_tim”, “gps_err”, “altitude”, “speed”

See also

get_flight_orientation_s2: get the quaternions of the flight path

Examples

Following the file and metadata structure of scihub:

>>> import os
>>> import numpy as np
>>> from s2d2.read_sentinel2 import list_central_wavelength_msi
>>> from s2d2.handler.sentinel2 import get_s2_image_locations

>>> s2_dir = '/data-dump/examples/'
>>> s2_name = 'S2A_MSIL1C_20200923T163311_N0209_R140_T15MXV_20200923T200821.SAFE'
>>> fname = os.path.join(s2_dir, s2_name, 'MTD_MSIL1C.xml')
>>> s2_df = list_central_wavelength_s2()

>>> s2_df, datastrip = get_s2_image_locations(fname, s2_df)
>>> path_det = os.path.join(s2_dir, s2_name, 'DATASTRIP', datastrip[17:-7])

>>> sat_tim, sat_xyz, sat_err, sat_uvw = get_flight_path_s2(path_det)

Notes

The metadata structure of MTD_DS looks like:

* MTD_DS.xml
└ n1:Level-1C_DataStrip_ID
   ├ n1:General_Info
   ├ n1:Image_Data_Info
   ├ n1:Satellite_Ancillary_Data_Info
   │  ├ Time_Correlation_Data_List
   │  ├ Ephemeris
   │  │  ├ GPS_Number_List
   │  │  ├ GPS_Points_List
   │  │  │  └ GPS_Point
   │  │  │     ├ POSITION_VALUES
   │  │  │     ├ POSITION_ERRORS
   │  │  │     ├ VELOCITY_VALUES
   │  │  │     ├ VELOCITY_ERRORS
   │  │  │     ├ GPS_TIME
   │  │  │     ├ NSM
   │  │  │     ├ QUALITY_INDEX
   │  │  │     ├ GDOP
   │  │  │     ├ PDOP
   │  │  │     ├ TDOP
   │  │  │     ├ NOF_SV
   │  │  │     └ TIME_ERROR
   │  │  └ AOCS_Ephemeris_List
   │  ├ Attitudes
   │  ├ Thermal_Data
   │  └ ANC_DATA_REF
   │
   ├ n1:Quality_Indicators_Info
   └ n1:Auxiliary_Data_Info

The following acronyms are used:

AOCS : attitude and orbit control system
CAMS : Copernicus atmosphere monitoring service
DEM : digital elevation model
GDOP : geometric dilution of precision
GPS : global positioning system
GRI : global reference image
IERS : international earth rotation and reference systems service
IMT : instrument measurement time
NSM : navigation solution method
NOF SV : number of space vehicles
PDOP : position dilution of precision
TDOP : time dilution of precision

basic reading of geo-referenced imagery

s2d2.mapping_input.read_geo_image(fname, boi=array([], dtype=float64), no_dat=None)[source]

This function takes as input the geotiff name and the path of the folder that the images are stored, reads the image and returns the data as an array

Parameters:

fname (string) – geotiff file name and path.
boi (numpy.array, size=(k,1)) – bands of interest, if a multispectral image is read, a selection can be specified
no_dat ({integer,float}) – no data value

Returns:

data ({numpy.array, numpy.masked.array}, size=(m,n), ndim=2) – data array of the band
spatialRef (string) – osr.SpatialReference in well known text
geoTransform (tuple, size=(8,)) – affine transformation coefficients.
targetprj (osgeo.osr.SpatialReference() object) – coordinate reference system (CRS)

See also

read_geo_info: basic function to get meta data of geographic imagery

Examples

>>> import os
>>> from s2d2.mapping_input import read_geo_image
>>> fpath = os.path.join(os.getcwd(), "data.jp2" )
>>> I, spatialRef, geoTransform, targetPrj = read_geo_image(fpath)

s2d2.mapping_input.read_geo_info(fname)[source]

This function takes as input the geotiff name and the path of the folder that the images are stored, reads the geographic information of the image

Parameters:

fname (string) – path and file name of a geotiff image

Returns:

spatialRef (string) – osr.SpatialReference in well known text
geoTransform (tuple, size=(8,1)) – affine transformation coefficients, but also giving the image dimensions
targetprj (osgeo.osr.SpatialReference() object) – coordinate reference system (CRS)
rows (integer, {x ∈ ℕ | x ≥ 0}) – number of rows in the image, that is its height
cols (integer, {x ∈ ℕ | x ≥ 0}) – number of collumns in the image, that is its width
bands (integer, {x ∈ ℕ | x ≥ 1}) – number of bands in the image, that is its depth

See also

read_geo_image: basic function to import geographic imagery data

coordinate transformations

s2d2.mapping_tools.ecef2llh(xyz)[source]

transform 3D cartesian Earth Centered Earth fixed coordinates, to spherical angles and height above the ellipsoid

Parameters:: xyz (np.array, size=(m,3), unit=meter) – np.array with 3D coordinates, in WGS84. In the following form: [[x, y, z], [x, y, z], … ]
Returns:: llh – np.array with angles and height. In the following form: [[lat, lon, height], [lat, lon, height], … ]
Return type:: np.array, size=(m,2), unit=(deg,deg,meter)

s2d2.mapping_tools.ecef2map(xyz, spatialRef)[source]

transform 3D cartesian Earth Centered Earth fixed coordinates, to map coordinates (that is 2D) in a projection frame

Parameters:

xyz (np.array, size=(m,3), float) – np.array with 3D coordinates, in WGS84. In the following form: [[x, y, z], [x, y, z], … ]
spatialRef (osgeo.osr.SpatialReference) – target projection

Returns:

xyz – np.array with planar coordinates, within a given projection frame

Return type:

np.array, size=(m,2), float

s2d2.mapping_tools.get_utm_zone(φ, λ)[source]

get the UTM zone for a specific location

Parameters:

ϕ_lim (float, unit=degrees) – latitude and longitude of a point of interest
λ_lim (float, unit=degrees) – latitude and longitude of a point of interest

Returns:

utm_zone – string specifying the UTM zone

Return type:

string

s2d2.mapping_tools.ll2map(ll, spatialRef)[source]

transforms angles to map coordinates (that is 2D) in a projection frame

Parameters:

llh (np.array, size=(m,2), unit=(deg,deg)) – np.array with spherical coordinates. In the following form: [[lat, lon], [lat, lon], … ]
spatialRef (osgeo.osr.SpatialReference) – target projection system

Returns:

xyz – np.array with 2D coordinates. In the following form: [[x, y], [x, y], … ]

Return type:

np.array, size=(m,2), unit=meter

Examples

Get the Universal Transverse Mercator (UTM) cooridnates from spherical coordinates:

>>> import numpy as np
>>> from osgeo import osr
>>> proj = osr.SpatialReference()
>>> proj.SetWellKnownGeogCS('WGS84')
>>> lat, lon = 52.09006426183974, 5.173794246145571# Utrecht University
>>> NH = True if lat>0 else False
>>> proj.SetUTM(32, True)

>>> xy = ll2map(np.array([[lat, lon]]), proj)
>>> xy
array([[ 237904.03625329, 5777964.65056734,       0.        ]])

s2d2.mapping_tools.map2ll(xy, spatialRef)[source]

transforms map coordinates (that is 2D) in a projection frame to angles

Parameters:

xy (np.array, size=(m,2), unit=meter) – np.array with 2D coordinates. In the following form: [[x, y], [x, y], … ]
spatialRef (osgeo.osr.SpatialReference) – target projection system

Returns:

ll – np.array with spherical coordinates. In the following form: [[lat, lon], [lat, lon], … ]

Return type:

np.array, size=(m,2), unit=(deg,deg)

transfer between mapping and image domain

s2d2.image_coordinate_tools.create_local_crs()[source]

create spatial refence of local horizontal datum

Returns:: crs – the coordinate reference system via GDAL SpatialReference description
Return type:: {osr.SpatialReference, string}

See also

s2d2.checking.mapping.is_crs_an_srs

s2d2.image_coordinate_tools.get_bbox(geoTransform, rows=None, cols=None)[source]

given array meta data, calculate the bounding box

Parameters:

geoTransform (tuple, size=(6,)) – georeference transform of an image.
rows (integer, {x ∈ ℕ | x ≥ 0}) – amount of rows in an image.
cols (integer, {x ∈ ℕ | x ≥ 0}) – amount of collumns in an image.

Returns:

bbox – bounding box, in the following order: min max X, min max Y

Return type:

np.ndarray, size=(4,), dtype=float

See also

map2pix, pix2map, pix_centers

Notes

Two different coordinate system are used here:

indexing   |           indexing    ^ y
system 'ij'|           system 'xy' |
           |                       |
           |       i               |       x
   --------+-------->      --------+-------->
           |                       |
           |                       |
image      | j         map         |
based      v           based       |

s2d2.image_coordinate_tools.map2pix(geoTransform, x, y)[source]

transform map coordinates to local image coordinates

Parameters:

geoTransform (tuple, size={(6,1), (8,1)}) – georeference transform of an image.
x (np.array, size=(m), ndim={1,2,3}, dtype=float) – horizontal map coordinate.
y (np.array, size=(m), ndim={1,2,3}, dtype=float) – vertical map coordinate.

Returns:

i (np.ndarray, ndim={1,2,3}, dtype=float) – row coordinate(s) in local image space
j (np.ndarray, ndim={1,2,3}, dtype=float) – column coordinate(s) in local image space

See also

pix2map

Notes

Two different coordinate system are used here:

indexing   |           indexing    ^ y
system 'ij'|           system 'xy' |
           |                       |
           |       i               |       x
   --------+-------->      --------+-------->
           |                       |
           |                       |
image      | j         map         |
based      v           based       |

s2d2.image_coordinate_tools.pix2map(geoTransform, i, j)[source]

transform local image coordinates to map coordinates

Parameters:

geoTransform (tuple, size={(6,1), (8,1)}) – georeference transform of an image.
i (np.array, ndim={1,2,3}, dtype=float) – row coordinate(s) in local image space
j (np.array, ndim={1,2,3}, dtype=float) – column coordinate(s) in local image space

Returns:

x (np.array, size=(m), ndim={1,2,3}, dtype=float) – horizontal map coordinate.
y (np.array, size=(m), ndim={1,2,3}, dtype=float) – vertical map coordinate.

See also

map2pix

Notes

Two different coordinate system are used here:

indexing   |           indexing    ^ y
system 'ij'|           system 'xy' |
           |                       |
           |       i               |       x
   --------+-------->      --------+-------->
           |                       |
           |                       |
image      | j         map         |
based      v           based       |

s2d2.image_coordinate_tools.pix_centers(geoTransform, rows=None, cols=None, make_grid=True)[source]

provide the pixel coordinate from the axis, or the whole grid

Parameters:

geoTransform (tuple, size={(6,), (8,)}) – georeference transform of an image.
rows (integer, {x ∈ ℕ | x ≥ 0}, default=None) – amount of rows in an image.
cols (integer, {x ∈ ℕ | x ≥ 0}, default=None) – amount of collumns in an image.
make_grid (bool, optional) – Should a grid be made. The default is True.

Returns:

X (np.ndarray, dtype=float) –
- “make_grid” : size=(m,n)
- otherwise : size=(1,n)
Y (np.ndarray, dtype=float) –
- “make_grid” : size=(m,n)
- otherwise : size=(m,1)

See also

map2pix, pix2map

Notes

Two different coordinate system are used here:

indexing   |           indexing    ^ y
system 'ij'|           system 'xy' |
           |                       |
           |       j               |       x
   --------+-------->      --------+-------->
           |                       |
           |                       |
image      | i         map         |
based      v           based       |

metric conversions

s2d2.unit_conversion.datenum2datetime(t)[source]

some write a date as a number “YYYYMMDD” conver this to numpy.datetime

Parameters:: t ({int,numpy.array}) – time, in format YYYYMMD
Returns:: dt – times in numpy time format
Return type:: {numpy.datetime64, numpy.array}, type=numpy.datetime64[D]

See also

datetime2datenum

s2d2.unit_conversion.datetime2calender(dt)[source]

Convert array of datetime64 to a calendar year, month, day.

Parameters:: dt (np.datetime64) – times of interest
Returns:: cal – calendar array with last axis representing year, month, day
Return type:: numpy array

s2d2.unit_conversion.datetime2datenum(dt)[source]

convert numpy.datetime to a date as a number “YYYYMMDD”

Parameters:: dt ({numpy.datetime64, numpy.array}, type=numpy.datetime64[D]) – times in numpy time format
Returns:: t – time, in format YYYYMMD
Return type:: {int, numpy.array}

See also

datenum2datetime

s2d2.unit_conversion.datetime2doy(dt)[source]

Convert array of datetime64 to a day of year.

Parameters:

dt (np.datetime64) – times of interest

Returns:

year (integer, {x ∈ ℕ}) – calender year
doy (integer, {x ∈ ℕ}) – day of year

s2d2.unit_conversion.deg2arg(θ)[source]

adjust angle to be in bounds of an argument angle

Parameters:: θ (unit=degrees) –
Returns:: θ
Return type:: unit=degrees, range=-180…+180

See also

deg2compass

Notes

The angle is declared in the following coordinate frame:

     ^ North & y
     |
- <--|--> +
     |
     +----> East & x

s2d2.unit_conversion.deg2compass(θ)[source]

adjust angle to be in bounds of a positive argument angle,like a compass

Parameters:: θ (unit=degrees) –
Returns:: θ
Return type:: unit=degrees, range=0…+360

See also

deg2arg

Notes

The angle is declared in the following coordinate frame:

     ^ North & y
     |
- <--|--> +
     |
     +----> East & x

s2d2.unit_conversion.deg2dms(ang)[source]

convert decimal degrees to degree minutes seconds format

Parameters:

ang ({float,np.array}, unit=decimal degrees) – angle(s) of interest

Returns:

deg ({integer,np.array}) – degrees
min ({integer,np.array}, range=0…60) – angular minutes
sec ({float,np.array}, range=0…60) – angular seconds

s2d2.unit_conversion.dms2deg(deg, min, sec)[source]

convert degree minutes seconds format to decimal degrees

Parameters:

deg ({integer,np.array}) – degrees
min ({integer,np.array}, range=0...60) – angular minutes
sec ({float,np.array}, range=0...60) – angular seconds

Returns:

ang – angle(s) of interest

Return type:

{float,np.array}, unit=decimal degrees

s2d2.unit_conversion.doy2dmy(doy, year)[source]

given day of years, calculate the day and month

Parameters:

year (integer, {x ∈ ℕ}) – calender year
doy (integer, {x ∈ ℕ}) – day of year

Return type:

day, month, year

handling formats and specifics

sentinel-2 metadata file structure

sentinel-2 tiling structure

s2d2.handler.mgrs.get_bbox_from_tile_code(tile_code, tile_path=None)[source]

Get the bounds of a certain MGRS tile

Parameters:

tile_code (string, e.g.: '05VMG') – MGRS tile coding
tile_path (string) – Path to the geometric metadata

Returns:

bounding box, in the following order: min max X, min max Y

Return type:

numpy.ndarray, size=(1,4), dtype=float

s2d2.handler.mgrs.get_geom_from_tile_code(tile_code, tile_path=None)[source]

Get the geometry of a certain MGRS tile

Parameters:

tile_code (string) – MGRS tile coding, e.g.: ‘05VMG’
tile_path (string) – Path to the geometric metadata

Returns:

Geometry of the MGRS tile, in lat/lon

Return type:

shapely.geometry.polygon.Polygon

See also

get_tile_code_from_geom, get_bbox_from_tile_code

s2d2.handler.mgrs.get_mgrs_tile(φ, λ)[source]

return a military grid reference system zone designation string. This zoning is used in the tiling of Sentinel-2.

Parameters:

ϕ (float, unit=degrees, range=-90...+90) – latitude
λ (float, unit=degrees) – longitude

Returns:

mgrs_code

Return type:

string

Notes

The MGRS tile structure is a follows “AABCC”

“AA” utm zone number, starting from the East, with steps of 8 degrees
“B” latitude zone, starting from the South, with steps of 6 degrees
“CC” 100 km square identifier

The following acronyms are used:

CRS : coordinate reference system
MGRS : US military grid reference system
WKT : well known text
UTM: universal transverse mercator projection

s2d2.handler.mgrs.get_tile_codes_from_geom(geom, tile_path=None)[source]

Get the codes of the MGRS tiles intersecting a given geometry

Parameters:

geom ({shapely.geometry, string, dict, GeoDataFrame, GeoSeries}) – geometry object with the given dict-like geojson geometry, GeoSeries, GeoDataFrame, shapely geometry or well known text, i.e.: ‘POLYGON ((x y, x y, x y))’
tile_path (string) – Path to the geometric metadata

Returns:

MGRS tile codes

Return type:

tuple

sentinel-2 metadata records

s2d2.handler.xml.get_array_from_xml(treeStruc)[source]: Arrays within a xml structure are given line per line Output is an array

generic

basic reading of Sentinel-2 data

advanced reading of Sentinel-2 sensors

basic reading of geo-referenced imagery

coordinate transformations

transfer between mapping and image domain

metric conversions

handling formats and specifics

sentinel-2 metadata file structure

sentinel-2 tiling structure

sentinel-2 metadata records

look at consistency

naming

mapping

data