Mapping

Fun NY Geographic Facts

|

The Town of Greenville is 39 sq miles, the Bourgh of Manhattan is 33 sq miles. A few more cows in Greenville though then Manhattan.

Staten Island about the size of Altona in Clinton County at 101 square miles. A few less wind turbines in Staten Island, and far less hot air.

Swimming Hole

All of New York City would fit into the area of Town of Long Lake, Webb, or Brookhaven.

The Town of Westerlo is approximately the same size as the Borough of Bronx at 58 square miles. Yes, but do they have Dick Rapp as their town supervisor?

The biggest town area-wise in New York State is Brookhaven on Long Island at 532 sq miles. Town of Webb (Old Forge) has 483 sq miles.

Barn Silloute Against TriMount

The difference between the year-round population of the Town of Webb and Town of Brookhaven is about the population of Albany County at 295,000 persons.

Yonkers is about the same size as Albany, NY at 21 square miles.

The typical town in Western NY is 36 square miles and is a nice square shape. Not all though as geography sometimes over rides that nice square shape.

Walk on By

The largest town in Western NY is Bath, NY at 101 square miles, or about the size of Staten Island.

Owego, near Binghamton is the largest Central New York town at 101 square miles. It should not be confused with Oswego which is near Syracuse.

People were asking how I was extracting the digital terrain model (tree/building elevation) from the NYS GIS LiDAR point clouds (LAS) files

|

People were asking how I was extracting the digital terrain model (tree/building elevation) from the NYS GIS LiDAR point clouds (LAS) files. While the digital surface models (bare earth) models are widely available in raster digital elevation model GeoTIFF, these files are bare earth and don’t include vegetation or buildings. But you can easily create raster digital elevation model GeoTIFF from the LIDAR point clouds with the PDAL – Point Data Abstraction Library.

Install PDAL

It’s easy if you are running Ubuntu Linux — it’s in standard distribution with 20.04 LTS or later. Most modern distributions of Linux include it in their repositories.

sudo apt install pdal

Download LiDAR Point Clouds

You can download the individual cloud files from NYS GIS FTP site. You can use the LiDAR Shapefile indexes on that website to figure out which file you need. Be aware you may have to download and store a significant amount of data — each roughly 556 acre tile (less then a square mile) is between 500 and 1,000 MB. So you will want to do this somewhere you have a lot of free unlimited data service and storage on your hard drive. And patience while it downloads.

Extract the Digital Terrain Model

You should create a text file called pdal_dtm.json or something similar. The text file should contain:

{ 
	"pipeline": [
		"18TWN220520.las",
		{ 
			"type": "filters.range",
			"limits": "Classification[1:1]"
		},
		{
			"filename":"dsm.tif",
			"gdaldriver":"GTiff",
			"resolution": 0.5,
			"output_type": "max",
			"type":"writers.gdal"
		}
	]
}

This creates the pdal “pipeline” configuration file for the conversion. An explanation of important lines you will need to change:

“18TWN220520.las” – Input point cloud LAS file
“type”: “filters.range” – Tells PDAL to filter points based on request
“limits”: “Classification[1:1]” – Layer to extract *
“filename”:”dsm.tif” – Export file name
“gdaldriver”:”GTiff” – Use GeoTIFF export (same as standard DEM files)
“resolution”: 0.5 – Resolution in meters for export, typically 0.5m or 1m depending on survey
“output_type”: “max” – Highest point reflected back within that point, to get digital surface elevation. You can also use min/max/mean/count/stdev/idw.
“type”:”writers.gdal” – Use the GDAL library for create GeoTIFF

* Layers in NYS GIS Point Clouds – “These point clouds will have at a minimum 2 classifications; Class 1 Unclassified, and Class 2 Bare Earth.” Layer 1 is the digital surface model, it contains building heights and tree/crop cover. Layer 2 “Bare Earth” is the same as what you get from the LiDAR Digital Elevation GeoTIFFs that are widely available on NYSGIS website. Be aware that the Class 2 points may have no-data areas, where no ground elevation was detected due to tree cover or other things blocking ground. The DEM files you download from NYSGIS have these areas filled in with interpolation from surrounding areas.

Once your configuration file is done, run the command:

pdal pipeline pdal_dtm.json

As the Point Cloud files are internally indexed, the export should take only a matter of seconds on a modern desktop computer. The exported DEM/GeoTIFF file will be between 15 to 25 MB, which is much smaller then the point cloud. Then you can load in Quantum GIS or your favorite non-free software GIS client like any other Raster DEM file.

To calculate the building or tree height, just subtract the Class 1 raster from Class 2 raster, using the Raster Calculator in QGIS or your favorite GIS program.  Or better yet, subtract Class 1 points from the already processed DEM files on NYSGIS website, as they have the missing bare earth data filled in.

How to deal with Shapefiles missing Projection information in QGIS

|

How to deal with Shapefiles missing Projection information in QGIS 🗺

Your best bet is to find the proper projection system by including the prj file in the same directory as the shapefile. If not available, try looking for the projection information (EPSG code) from where you got the file or the included XML meta data which can be opened in your web browser.

As a last resort, you can use QGIS Find Projection tool (use the Type to locate box in the lower left) after loading the shapefile to get potential projections. That won’t get you one answer but will narrow down the options. Be aware that using the wrong datum with the right coordinate system will introduce small errors in the projection – often a couple hundred meters or less, so if things are slightly off try a different datum. But try to find the official projection from the prj or xml file.

Store your own data as geopackages to avoid lost prj files and future projection issues.

How to label a point in QGIS based on a vector layer below it

|

How to label a point in QGIS based on a vector layer below it

For example, you might have a point layer and want to label what town it is in. You might have the Census Tiger County Subdivision layer named ‘tl_2020_36_cousub20’, with the field ‘NAMELSAD20’ that you want to return from your point or polygon layer.

attribute(array_first(overlay_nearest(‘tl_2020_36_cousub20′,$currentfeature)),’NAMELSAD20’)

This will give you the nearest feature.

You need to wrap it in array_first as shown above, as overlay_nearest returns an array with one element. Then attribute just returns the attribute based on the object provided.

Alternatives include overlay_contains, overlay_crosses, overlay_intersects, overlay_touches, overlay_within which will return what you want. Try searching the QGIS Help in the QGIS Expressions builder.

You can also use this style of expression to create a spatial join or nearest neighbor join, inside the field calculator.

I’ve noticed that QGIS doesn’t always do a good job simplifying geometries when exporting to an SVG file

|

I’ve noticed that QGIS doesn’t always do a good job simplifying geometries when exporting to an SVG file. While you can manually simplify layers, an easier method is as follows:

1) Open Layer Styling

2) Go to Symbology

3) Switch fill to Geometry Generator

4) Use the Simplify Expression to Simplify on the Fly, such as simplifying by 100 map units, which when projected in UTM would be meters.

simplify($geometry,100)