How to get point-in-polygon for different time stamps?

Question

I have an input polygon vector layers which is essentially a single KML file containing different areas over Europe. I also have an input point vector layer which is actually flight data, which contains Lat/Lon and time stamp (24hr worth of data) of aircrafts.

What I am trying to do is to get the maximum number of aircrafts that will be inside each of these polygons over the 24hr window.

I have tried Points in Polygon, but that gives me the total number of points registered inside each polygon, but they may occur at different times of the day. What I basically want is a max of max.

EDIT

@Henhuy when i include: print(max_points, max_id, feature.id()) in the loop, i get the max values printed as follows:

But this is how it is shown in the attribute table:

In my case, row 884 (object id 885) should have a max_point of 4, but it is showing this as NULL.

Answer 1

I tried it successfully with following code. It counts maximum points per polygon per day. For my code, point layer needs attribute "date"; then, for each date occurrences are count in dictionary and stored in layer for every polygon at end of loop. Note: New column is added to polygon layer called 'max_points' and filled with results

Code is this (I added it to python script editor in QGis, loaded both layers via selecting each layer and calling iface.activeLayer() then called the function with both layers:

Updated Code (to store values in layer):

def max_points_in_polygon_per_day(poly_layer, points_layer):
    from PyQt4.QtCore import QVariant
    poly_layer.startEditing()
    newField = QgsField('max_points', QVariant.Int)
    poly_layer.dataProvider().addAttributes([newField])
    poly_layer.updateFields()
    max_id = poly_layer.getFeatures().next().fieldNameIndex('max_points')
    for feature in poly_layer.getFeatures():
        geom = feature.geometry()
        days = {}
        for point_feature in points_layer.getFeatures():
            point_geom = point_feature.geometry()
            if point_geom is None:
                continue
            if point_geom.within(geom):
                date = point_feature.attribute('date')
                if date in days:
                    days[date] += 1
                else:
                    days[date] = 1
        if len(days) == 0:
            max_points = 0
        else:
            max_points = max(days.values())
        poly_layer.changeAttributeValue(feature.id(), max_id, max_points)
    poly_layer.commitChanges()

Looping through both feature lists could take a will and maybe is not the fastest solution - feel free to optimize my code!

EDIT

Here a short explanation how to get it to work.

First click on Python-Console:

Then open python editor:

Insert code there and click on "Run script" to reload script. Afterwards click on point layer and type following into python console:

point_layer = iface.activeLayer()

Then click on polygon layer and type:

poly_layer = iface.activeLayer()

Then type:

max_points_in_polygon_per_day(poly_layer, point_layer)

And it should print out something useful...

Answer 2

Using the plugin 'Points2One' or 'PointsToPath' you can create lines from your points by grouping the features together on the aircraft registration number. The end product should be each flight path.

Next intersect the lines over the polygons and finally extract the centroids from the intersected line segments. Using 'Count Points in Polygon' will give you the total number of aircraft within each polygon.