Debugging drill: clean a messy application

A fictional application arrives from its applicant: two parcels in one shapefile, with every problem of the area data checklist once, the problems area data from the field really has. You find each one with the checklist, fix it with QGIS’s own tools, and check the fix: each check shows the problem before its fix and passes after it. At the end, the application’s corners are its technical description’s corners again, to the millimeter.

ONE WAY, NOT THE ONLY WAY

The fixes here are examples, tested with the drill’s data. Other tools are equally valid, as long as the checks pass. Where a choice affects correctness (reprojecting a layer instead of assigning it a CRS), the drill says why.

Before we begin

NOTE

For this drill, you will use the d1-messy-application subfolder, and bnhr_qgis-sflma-mapping-001.gpkg, in the data subfolder of the BNHR_QGIS-for-SFLMA-Mapping folder.

What we want to achieve

1. Find every problem of the application’s area data with the area data checklist.
2. Fix each one, in an order where each fix makes the next one possible.
3. Show that the cleaned application is its technical description: the same corners, in number and position, and the same areas.

What we have

1. The application of the fictional Piñahan sa Lambak Agroforestry Cooperative in General Nakar, Quezon, as delivered, in d1-messy-application:
   • d1_application.shp (with its .shx, .dbf and .prj) - the application area as the applicant sent it
   • d1_td.csv and d1_td_corners.csv - its technical description, in PRS 1992 Philippines Zone III: the lines and the corners of each parcel
   • README.txt - its parcels: parcel 1, 7 corners and 14.2783 hectares; parcel 2, 5 corners and 5.4771 hectares
2. ref_barangays and ref_municipalities in bnhr_qgis-sflma-mapping-001.gpkg
3. A project in EPSG:3123 - PRS92 / Philippines zone 3, as set up in Exercise 1.1

Relevant QGIS knowledge/skills

• The area data checklist
• Validating the geometry
• Computing the area and location attributes
• The QGIS Processing Framework

Finding the problems

D1.1. Inspecting the application

Level: Beginner

1. Load d1_application.shp, ref_barangays, and ref_municipalities.
2. Go through the checklist, and write down what each check reads. Don’t fix anything yet.

Two parcels, one of them across the boundary with Infanta, and four small holes

Checks 10 and 12 are measured in meters, so read them after D1.3, once the layer is in the national grid. What you should have found:

#	Check	The application reads	Fixed in
1	Where it came from	249 vertices for the 12 corners of its technical description, and holes the size of raster cells: traced, not plotted	D1.6, D1.8, shown in D1.12
2	Files	No .cpg	D1.2, D1.12
3	CRS	EPSG:4326 - WGS 84; the technical description is in PRS 1992 Philippines Zone III	D1.3
4	Encoding	Pi�ahan sa Lambak Agroforestry Cooperative	D1.2
5	Field names	MUNICIPALI, cut at 10 characters; Brgy, not the guide’s name	D1.9
6	Geometry type	Polygon (MultiPolygonZ), every Z value 0	D1.4
7	Parcels	One feature of two parcels	D1.5
8	Validity	Invalid: Ring self-intersection	D1.7
9	Holes	Four: three in parcel 1, one in parcel 2	D1.6
10	Corners	137 and 87 vertices for 7 and 5 corners; a line of 0.40 m	D1.8
11	Attributes	REGION CALABARZON; MUNI_CTY General Nacar; AREA the whole application’s	D1.10
12	Location	0.6139 ha of parcel 2 in Infanta	D1.11, flagged

The fixes come in this order: first what lets you read the data (its names, then its coordinates in meters), then the geometry from the whole to its details (Z values, parts, holes, vertices), then the fields and their values, which the geometry’s changes make out of date, and last the check against the technical description.

Fixing the problems

D1.2. Setting the encoding

Level: Beginner

1. Open the attribute table of d1_application. The application’s name reads Pi�ahan: the shapefile has no .cpg, so QGIS reads its text as UTF-8, and the program that wrote it used Windows’ encoding.

The ñ read as �, a field cut at 10 characters, and a field not named as the guide names it

2. Open Layer Properties ‣ Source and set Data source encoding to windows-1252, then click OK.

Data source encoding: windows-1252

Check: the attribute table now reads Piñahan sa Lambak Agroforestry Cooperative. The setting holds only in this project: D1.12 saves the cleaned layer in UTF-8, with a .cpg, so that it reads correctly everywhere. Every tool you run on the layer from now on reads the names correctly.

D1.3. Reprojecting to the technical description’s zone

Level: Beginner

1. Open Layer Properties ‣ Information. The layer is in EPSG:4326 - WGS 84, its extent in degrees, while its technical description is in PRS 1992 Philippines Zone III.

Encoding UTF-8 (no .cpg), Polygon (MultiPolygonZ), and WGS 84

2. Run Processing Toolbox ‣ Reproject layer on d1_application with the Target CRS EPSG:3123 - PRS92 / Philippines zone 3. QGIS converts the coordinates with the transformation PRS92 to WGS 84 (1).

Reproject layer to EPSG:3123

REPROJECT, DON'T ASSIGN

The Assigned Coordinate Reference System on the Source page (or Set CRS) only relabels the coordinates, as its own note says. Assigning EPSG:3123 to this layer reads its degrees as meters, and puts the application 1,740 km from General Nakar.

Check: the Reprojected layer is in EPSG:3123, and lies on ref_barangays in General Nakar, beside Infanta.

D1.4. Dropping the Z values

Level: Beginner

1. The Information page read Polygon (MultiPolygonZ): the layer is 3D. Preview this expression in the Field Calculator () on Reprojected: it gives 0, so the Z values carry nothing, and nothing is lost by dropping them.

z_max($geometry)

2. Run Processing Toolbox ‣ Drop M/Z values on Reprojected, with Drop M Values and Drop Z Values checked.

Check: the output’s Layer Properties ‣ Information reads Polygon (MultiPolygon).

D1.5. One feature per parcel

Level: Beginner

1. The layer has one feature, and this expression gives 2 for it: its two parcels are one feature.

num_geometries($geometry)

2. Run Processing Toolbox ‣ Multipart to singleparts on the output of D1.4.

Check: the output, Single parts, has two features, and the expression gives 1 for each. Both have the whole application’s attributes, AREA among them: D1.10 computes them again.

D1.6. Deleting the holes

Level: Intermediate

1. This expression gives 3 for parcel 1 and 1 for parcel 2:

num_interior_rings($geometry)

2. Zoom to the holes. Each is a square of 10 m, 100 m², except one of two squares touching at a corner, 200 m²: the cells of a raster the parcels were traced from. A technical description describes the boundary only, and this one has no exclusions, so the holes are artifacts, not exclusions.
3. Write down each parcel’s area, from the preview of area($geometry).
4. Run Processing Toolbox ‣ Delete holes on Single parts, with Remove holes with area less than 500 square meters: above the largest artifact, so that a larger hole, which could be a real exclusion, would stay.

Remove holes with area less than 500 m²

Check: the expression gives 0 for both parcels of Cleaned. Parcel 1 grew by 400 m² and parcel 2 by 100 m², 0.0500 ha in all: record it, for the data’s owner.

ASK BEFORE YOU DELETE A HOLE

A hole that is a real exclusion (a road, a settlement, another tenure) is not an error, and deleting it applies for land the applicant did not ask for. Decide with the data’s owner which holes are artifacts.

D1.7. Checking the validity again

Level: Beginner

1. In D1.1, Processing Toolbox ‣ Check validity (GEOS) found parcel 1 invalid, with a Ring self-intersection: at the corner where the two cells of one hole touch, the hole’s ring passes through the same point twice.

Ring self-intersection, where the two cells of one hole touch

2. The error was in a hole, and D1.6 deleted the hole: this problem has no fix of its own. Run Check validity again on Cleaned.

Check: both parcels are valid.

DON'T FIX IT WITH FIX GEOMETRIES

Processing Toolbox ‣ Fix geometries also makes parcel 1 valid, by splitting that hole in two: parcel 1 then has four holes, and every hole stays. As in Exercise 2.5, fix an error where it comes from.

D1.8. Keeping only the corners

Level: Intermediate

1. On Cleaned, preview the number of vertices of each parcel: 137 for parcel 1 and 87 for parcel 2, for the technical description’s 7 and 5 corners.

num_points($geometry) - 1

2. Preview the shortest line between two vertices: 0.40 m.

array_min(array_foreach(geometries_to_array(segments_to_lines($geometry)), length(@element)))

A vertex every 10 m on each line of parcel 1

3. Every extra vertex is on a straight line between two corners, the 0.40 m one included. Run Processing Toolbox ‣ Simplify on Cleaned, with Simplification method Distance (Douglas-Peucker) and Tolerance 0.01 meters.

Simplify removes only the vertices within 1 cm of a straight line

Check: on Simplified, the vertices are 7 and 5, the shortest line is 161.02 m, and area($geometry) is what it was.

Parcel 1’s seven corners, and nothing else

D1.9. Giving the fields the guide’s names

Level: Intermediate

1. The attribute table has MUNICIPALI, the applicant’s MUNICIPALITY cut at 10 characters, and Brgy, which is not the guide’s name: the guide’s fields are MUNI_CTY and BARANGAY (Exercise 2.1).
2. Run Processing Toolbox ‣ Refactor fields on Simplified. Make each row one of the guide’s fields, in the guide’s order, with its type and length: the Source Expression of MUNI_CTY is "MUNICIPALI", and of BARANGAY "Brgy"; the others keep their own field.

MUNI_CTY from MUNICIPALI and BARANGAY from Brgy, the guide’s lengths

3. Name the output d1_clean in the Layers panel.

PRO-TIP

If you kept the empty shapefiles of Exercise 2.1 as templates, choose the area boundary in Load fields from template layer and click Load Fields: Refactor fields then lists the guide’s fields, and you only set the sources of the two that differ.

Check: d1_clean has the fields AREA, NAME, TENURE, REGION, PROVINCE, MUNI_CTY, and BARANGAY, in that order.

D1.10. Fixing the attributes

Level: Intermediate

1. Three values are wrong: REGION CALABARZON, which Exercise 5.4’s REGION expression selects () on both parcels; MUNI_CTY General Nacar, spelled otherwise than ref_barangays spells it; and AREA 19.7054, the whole application’s as delivered, on both parcels.
2. With the Field Calculator, update REGION to 'IV-A', and AREA, MUNI_CTY, and BARANGAY with the expressions of Exercise 2.6.

MUNI_CTY spelled as the reference layer spells it

Check: the REGION expression selects nothing, and the attributes read:

Each parcel’s area, the region in Roman numerals, and the municipality spelled right

D1.11. Checking the location

Level: Intermediate

1. Preview the hectares of each parcel outside General Nakar on d1_clean: 0 for parcel 1, and 0.6139 for parcel 2.

area(difference($geometry, aggregate('ref_municipalities', 'collect', $geometry,
filter:="muni_city" = 'General Nakar'))) / 10000

2. To see where, select General Nakar in ref_municipalities and run Processing Toolbox ‣ Difference with d1_clean as the input, ref_municipalities as the overlay, and Selected features only checked for the overlay.

Parcel 2 crosses into Infanta, by 0.6139 ha on the reference layer

This is not yours to fix. The parcels are the technical description’s (D1.12 shows it), so either the application does reach into Infanta, and its MUNI_CTY should list both municipalities, or the reference layer’s boundary is off there: reference layers are indicative. MUNI_CTY reads only General Nakar because ref_barangays has General Nakar’s barangays only.

Check: your findings say, for the data’s owner: parcel 2 has 0.6139 ha in Infanta, on GADM’s boundary.

Checking the result

D1.12. Checking against the technical description, and saving

Level: Intermediate

1. Load d1_td_corners.csv as points, as in Exercise 2.3 (X field EASTINGS, Y field NORTHINGS, EPSG:3123).
2. On d1_td_corners, Select by Expression the corners that are not on a vertex of d1_clean. No corner should be selected.

distance($geometry, nodes_to_points(aggregate('d1_clean', 'collect', $geometry))) > 0.001

The technical description’s corners that are more than 1 mm from a vertex of the cleaned application

Every corner on a vertex, and every vertex a corner

3. Right-click d1_clean ‣ Export ‣ Save Features As…: format ESRI Shapefile, file name d1_clean_application.shp, encoding UTF-8.

Check: as delivered, the application had 249 vertices for the technical description’s 12 corners. Now no corner is selected, the parcels have 7 and 5 vertices, and their AREA is the README’s, 14.2783 and 5.4771 hectares. d1_clean_application has its .cpg, and its name reads correctly with no encoding set.

Stretch: the cleaning as a model

D1.13. Saving the cleaning as a Processing model

Level: Advanced

An office cleans many applications. A Processing model saves the cleaning once, as one tool, so that the next application is cleaned in one run, with the same settings, and the model records what was done to the data. Read about models in QGIS: Essentials - The QGIS Processing Framework.

1. Open Processing ‣ Graphical Modeler, and name the model Clean an application.
2. Add its inputs, from the Inputs panel: a Vector Layer (polygon) named Application area; a CRS named Target CRS, by default EPSG:3123; and two Numbers, Remove holes smaller than (m²), 500, and Simplify tolerance (m), 0.01.
3. Add the algorithms of D1.3 to D1.10, from the Algorithms panel, each taking the previous one’s output: Reproject layer (its target CRS the Target CRS input), Drop M/Z values (both checked), Multipart to singleparts, Delete holes (its area the input), Simplify (Douglas-Peucker, its tolerance the input), and Refactor fields, as in D1.9 but with the values of D1.10 as its source expressions: 'IV-A' for REGION, and the expressions of Exercise 2.6 for AREA, MUNI_CTY, and BARANGAY. Name Refactor fields’ output Clean application.

Six steps in a column, each input beside the step it sets

4. Save the model, set d1_application’s encoding as in D1.2, and run the model on it.

Check: its output passes the checks of D1.12.

The model runs every step every time, and a step changes nothing where nothing is wrong: reprojecting a layer already in the target CRS, or simplifying a parcel with only its corners, leaves it as it is. The choices are its inputs. What it can’t do is judge: whether holes are artifacts (its hole area keeps any larger hole), the encoding to read the input in, and the location, which stays a finding for the data’s owner. Its Refactor fields is this application’s: change its sources for an applicant who names the fields otherwise.