0

I am trying to create a transect profile and extract information for each specific bands of Landsat which will be further used to create a transect graph. Is there any direct way of achieving this using QGIS?

1

Is this a One time Task?

YES:

You should look on "points along geometry" tool followed by "raster sampling" and further "points to path" all available on the processing tools.

No, it is a step on a much larger workflow.

Then use python (script or console). There you can iterate over the points on a linestring feature, and get the raster value on each point. Or you can make regular sampling on stright lines, or you can... anything you want.

I once did this with something like the following, (and worked):

dep = os.path.join(DIR,FILE) # get file path
dep_lyr = QgsRasterLayer(dep,'some-name')

uri = 'MultiPointZ?{}'.format('&'.join(['crs=epsg:4326','field=linea:integer',
                                        ,'index=yes'])) #you may have more fields.
grilla = QgsVectorLayer(uri,'muestreo','memory') # Layer to store results
campos = grilla.dataProvider().fields()

for feat in perfiles.getFeatures():
    # perfiles was an ad-hoc linestring layer.
    raster = dep_lyr
    xmax = raster.extent().xMaximum()
    dp = raster.dataProvider()
    #
    ndivs = int(length / sample_distance) #had precalculated values
    #
    geom = QgsGeometry(feat.geometry())
    # check if you need to transform your features coordinates, 
    # When using the python console there are no automatic reprojections
    ml = geom.get() #I had multilinestring features, 
                    # you can have linestrings instead.
    with edit(grilla):
        for nl,linea in enumerate(ml):
            # Iterate over linestrings in multilinestring.
            mp = QgsMultiPoint()
            subfeat = QgsFeature(campos)
            subfeat['linea'] = nl # set line number inside linestring
            # Regular sample on stright lines.
            X = np.linspace(linea.xAt(0),linea.xAt(1),ndivs)
            Y = np.linspace(linea.yAt(0),linea.yAt(1),ndivs)
            for x,y in zip(X,Y):
                _ = mp.addGeometry(QgsPoint(x,y,0))

            for pt in mp:
                #iterate over points 
                if xmax > 180:
                # do some dateline wrapping, you probably don't need it
                    if pt.x() < 0:
                        s = QgsPointXY(pt.x()+360,pt.y())
                    else:
                        s = QgsPointXY(pt)
                else:
                    s = QgsPointXY(pt)

                # Here is where we sample the raster.
                _ = pt.setZ(dp.sample(s,1)[0])
                #        
            # add geometry to feature
            _ = subfeat.setGeometry(mp)
            # add feature to layer
            _ = grilla.addFeature(subfeat)

Now you have one MultiPointZ feature for each profile. You should probably Read the docs on Raster Sampling for the band selection step.

Exporting the layer to any text format will give you a nice wkt 3d point list.

Or you can export your data from the console to a text file instead of creating a new layer.

| improve this answer | |
1

You may create transect profile using 'Profile Tool' that can be directly installed if you are connected with internet. You can install this plugin through 'Manage and Install plugins option available in 'Plugins' menu.

Transect

Resultant Profile for six bands

| improve this answer | |

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.