5

I'm trying to find either a list of all of the MODIS Terra swath footprints, or boundary dataset (i.e. polygons) with the footprint code as an attribute. I've searched high and low on various NASA and MODIS-related sites, but I can't find that specific data anywhere.

Does anyone know where I could find this information?

1

You can find a shape file with this information under the following link:

http://book.ecosens.org/modis-sinusoidal-grid-download/

It contains the horizontal and vertical tile numbers as attributes 'h' and 'v'.

| improve this answer | |
  • 1
    Sorry, but the MODIS Sinusoidal tiling system does not apply to Swath products, which is what @tinyplanet00 explicitly asked for. @Mr.Che's answer taking the official /geoMeta/ datasets into account is more target-oriented. – fdetsch Sep 18 '18 at 12:42
2

Also you could try to download txt files with footprint corners coordinates located here: https://ladsweb.modaps.eosdis.nasa.gov/archive/geoMeta/6/TERRA/

File example: https://ladsweb.modaps.eosdis.nasa.gov/archive/geoMeta/6/TERRA/2018/MOD03_2018-01-01.txt

Python 3 Code to convert *.txt MODIS footprints to ESRI shapefile (as @mikoontz requested). QGIS 3 required.

Disadvantage: the problem of the 180-degree meridian is not solved - if a polygon lies on it, the geometry becomes incorrect (it looks like a huge footprint crosses the Earth).

enter image description here

import qgis.core
import qgis.PyQt.QtCore


def parse_txt_to_dict(txt_filepath):
    txt_as_dict = []
    with open(txt_filepath) as file:
        string = file.read()

    string = string.replace('\\n', '\n')
    lines_list = string.split('\n')

    # catch filed names
    fields_list = None
    for line in lines_list:
        if line.startswith("#") and 'GranuleID' in line:
            line = line.replace('# ', '')
            fields_list = line.split(',')
        # print(len(fields_list), fields_list)

    # catch values
    if fields_list is None:
        print('fields_list is None')
        return None

    for line in lines_list:
        if not line.startswith("#"):
            line = line.split(',')
            if not line == ['']:
                dictionary = dict(zip(fields_list, line))
                txt_as_dict.append(dictionary)

    return txt_as_dict


def dict_to_shp(dicts_list, shp_path):
    if dicts_list == []:
        print('Empty list. Skip.')
        return None

    # fields
    keys = dicts_list[0].keys()
    fields = qgis.core.QgsFields()
    for key in keys:
        fields.append(qgis.core.QgsField(key, qgis.PyQt.QtCore.QVariant.String))

    layer_name = 'test'
    # layer = qgis.core.QgsVectorLayer("Polygon?crs=epsg:4326", layer_name, "memory")
    qgis.core.QgsVectorLayer("Polygon?crs=epsg:4326", layer_name, "memory")

    writer = qgis.core.QgsVectorFileWriter(vectorFileName=shp_path,
                                           fileEncoding='UTF-8',
                                           fields=fields,
                                           geometryType=qgis.core.QgsWkbTypes().Polygon,
                                           srs=qgis.core.QgsCoordinateReferenceSystem("EPSG:4326"),
                                           driverName='ESRI Shapefile',
                                           )
    # geometry
    for dictionary in dicts_list:
        points = [qgis.core.QgsPointXY(float(dictionary['GRingLongitude1']), float(dictionary['GRingLatitude1'])),
                  qgis.core.QgsPointXY(float(dictionary['GRingLongitude2']), float(dictionary['GRingLatitude2'])),
                  qgis.core.QgsPointXY(float(dictionary['GRingLongitude3']), float(dictionary['GRingLatitude3'])),
                  qgis.core.QgsPointXY(float(dictionary['GRingLongitude4']), float(dictionary['GRingLatitude4'])),
                  qgis.core.QgsPointXY(float(dictionary['GRingLongitude1']), float(dictionary['GRingLatitude1'])),
                  ]

        feature = qgis.core.QgsFeature()
        feature.setGeometry(qgis.core.QgsGeometry.fromPolygonXY([points]))
        attr_vavalues = list(dictionary.values())
        feature.setAttributes(attr_vavalues)

        writer.addFeature(feature)
    del writer


def main():
    txt_path = r'C:\DELETE\GIS SE\txt\MOD03_2018-01-01.txt'
    shp_path = r'C:\DELETE\GIS SE\shp\MOD03_2018-01-01_shp.shp'
    dicts_list = parse_txt_to_dict(txt_path)
    dict_to_shp(dicts_list, shp_path=shp_path)


if __name__ == "__main__":
    main()
    print('Done!')
| improve this answer | |
  • This should be the correct answer, I think. The README for the geoMeta product says it's not an official MODIS product (though distributed by the same folks), so it isn't searchable! Nice find! I'd be interested in seeing your approach to converting the corner boundaries to polygons if you ever happen to read this! – mikoontz Feb 15 at 21:31
  • @mikoontz Answer updated. – Comrade Che Feb 18 at 5:26

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.