FlightGear wiki - User contributions [en]

User:Wlbragg

2026-04-01T23:14:13Z

Wlbragg: /* Mission Generator Addon */

[[File:Kansas Scenery Three.jpg|632x632px|right|ElDorado, Kansas - ElDorado Lake]]
== Kansas/Midwest Scenery Development ==
<big>Hello and welcome to my Kansas/Midwest scenery development page.</big>

Hopefully in the weeks, months and years to come, I can contribute some useful Midwest Scenery and Airport Data.

=== Airport Layouts ===
So far, I have updated and submitted airport layouts to the [[FlightGear]] project for
*KEQA
*3AU
*KRCP
I have a few more in the works and plan to eventually add terminals and structures to as many Kansas Airports as I have the motivation to finish.

[[File:KEQA_Airport.jpg|350px|KEQA]]

=== Kansas Custom Scenery ===
I also have custom scenery generated from a combination of
*NLCD (land cover)
*OSM (roads, rivers and rail)
*SRTM (water data)
*New CS Data
*custom textures
for the entire state of [https://1drv.ms/u/s!AmL-8HrJKQ8ciRZl-sKR0HWnKPJQ?e=YuJoQ4 Kansas]

This is currently very heavy data but works surprisingly well on my mid level CPU and graphics card.

<gallery caption="Kansas Scenery" mode="packed-hover">
Kansas Scenery One.jpg|ElDorado, Kansas
Kansas Scenery Four.jpg|ElDorado Lake
Kansas Scenery Two.jpg|Kansas City
KICT.jpg|KICT - Wichita Dwight D. Eisenhower National Airport (Mid-Continent)
KIAB.jpg|KIAB - McConnell AFB, Wichita
Kansas.jpg|Kansas Sunset
KMCI.jpg|KMCI - Kansas City International Airport
ElDorado Lake Sunset.jpg|Sunset at ElDorado Lake
Cub Over Eldorado Lake.jpg|Cub at ElDorado Lake
</gallery>

=== Ohio Custom Scenery ===
Custom scenery generated from a combination of
*CS Data based on NLCD (land cover)
*OSM (roads, rivers and rail)
*SRTM (water data to come)
*custom textures
for the entire state of [https://1drv.ms/u/s!AmL-8HrJKQ8ciAqGzBvP5O-kXZRz?e=1i4xWm Ohio].

<gallery caption="Ohio Scenery" mode="packed-hover">
Cincinnati.jpg|KLUK - Cincinnati Municipal Airport
Toledo-Harbor View.jpg|Toledo - Harbor View
Lake_Erie_2.jpg|KCLE Cleveland-Hopkins International Airport
KPCW.jpg|KPCW - Carl R Keller Field
South Lorain.jpg|KCAK Akron-Canton Regional Airport
KBKL.jpg|KSKY - Griffing Sandusky Airport
Lake Erie.jpg|KCMH Port Columbus International Airport
KDAY_James_M._Cox_Dayton_International_Airport.jpg|KDAY James M. Cox Dayton International Airport
KYNG Youngstown-Warren Regional Airport.jpg|KYNG Youngstown - Warren Regional Airport
Near Youngstown.jpg|Near Youngstown
</gallery>

== Addon Development ==
[[File:Cargo-position.jpg|thumb|501x501px|Cargo placement UI.]]

=== Cargo Towing Addon ===
Cargo Towing is now in the form of an addon. The Cargo Towing addon adds the ability for any helicopter to haul cargo from either a slung position. The AirCrane retains a docked cargo capability. Now included, hooks for adding load weight and winched rope length to any aircraft that has the instrumentation do display this data.. Using the Wildfire addon in addition to Cargo Towing addon, you can use a bambi bucket for firefighting. Also added a special feature that allows the AirCrane to load water from any of the 3 Water Tank cargos. The addon extended the sling load Cargo Towing capability to any aircraft.

[https://svn.code.sf.net/p/flightgear/fgaddon/trunk/Addons/CargoTowingAddon Cargo Towing Addon]

[[File:Intercept1.jpg|left|thumb|357x357px]]

=== Mission Generator Addon ===
The Mission Generator addon is a continuing project that is designed to give FlightGear some ability to easily generate different mission orientated challenges. So far there are currently ...

* Wildfire events. Using Wildfire and Cargo Towing addons you have the capability to put these fire out using various aircraft and methods. Including Water drop with tankers, bambi bucket or the Aircrane's jet spray boom apparatus.
* Search and rescue, water and land, including life flight to hospital.
* Intercept missions, including an airspace intrusion

All missions are fully customizable and can be either on demand or set a random range and time for the event, per event type.

[https://svn.code.sf.net/p/flightgear/fgaddon/trunk/Addons/MissionGenerator Mission Generator Addon]

[[File:Intercept3.jpg|thumb|400x400px|Vectored to target.|left]]
[[File:Intercept2.jpg|left|thumb|408x408px|Target data vectoring.]][[File:Misson Generator Help.png|thumb|Help visual for the Mission Generator UI|border|center|820x820px]]
=== WildfireAddon ===
[[File:Bambi-bucket-full.jpg|thumb|761x761px|Bucket fills by submerging in any water source.]]The Wildfire addon originally was a FlightGear stand alone nasal module written and developed by Anders Gidenstam and incorporated in fgdata

Credit:

Wildfire - A cellular automaton forest fire model with the ability to spread over the multiplayer network.Written and developed by Anders Gidenstam (anders(at)gidenstam.org) Copyright (C) 2007 - 2012 Anders Gidenstam (anders(at)gidenstam.org) Modifications and conversion to Addon by Wayne Bragg wlbragg<blockquote>

The cellular automata model used here is loosely based on A. Hernandez Encinas, L. Hernandez Encinas, S. Hoya White, A. Martin del Rey, G. Rodriguez Sanchez, "Simulation of forest fire fronts using cellular automata", Advances in Engineering Software 38 (2007), pp. 372-378, Elsevier.</blockquote>[https://svn.code.sf.net/p/flightgear/fgaddon/trunk/Addons/Wildfire Wildfire Addon]
[[File:Bambi-bucket-dump.jpg|thumb|455x455px|Dumping the bucket will suppress the wildfire..|alt=|none]]

[[File:Bambi-busket2.jpg|thumb|741x741px|View from engineering station.|alt=|none]]

== Future Endeavors ==

=== Surprise ===

[[File:Iss4.jpg|512px|International Space Stay|alt=|left|thumb]][[File:Aircrane Over the Mojave.jpg|1024px|AirCrane over the Mojave Desert in the Palmdale/Landcaster area. Compositor shadows, OSM2CITY models, ALS lighting and Orthographics texture.]]

[[File:Aircrane Over the Mojave 2.jpg|1024px|AirCrane over the Mojave Desert in the Palmdale/Landcaster area. Compositor shadows, OSM2CITY models, ALS lighting and Orthographics texture.]]
{{-}}

Watch this page for future releases of various project for [[FlightGear]]. I hope to have a couple of surprises just around the bend.

File:Misson Generator Help.png

2026-04-01T23:11:08Z

Wlbragg: Wlbragg uploaded a new version of File:Misson Generator Help.png

=={{int:filedesc}}==
{{Information
|description={{en|1=Help visual for the Mission Generator UI}}
|date=2026-04-01
|source={{own}}
|author=[[User:Wlbragg|Wlbragg]]
|permission=
|other versions=
}}

=={{int:license-header}}==
{{self|cc-by-sa-4.0}}

User:Wlbragg

2026-04-01T19:48:46Z

Wlbragg: Added a help image for the Mission Generator UI

[[File:Kansas Scenery Three.jpg|632x632px|right|ElDorado, Kansas - ElDorado Lake]]
== Kansas/Midwest Scenery Development ==
<big>Hello and welcome to my Kansas/Midwest scenery development page.</big>

Hopefully in the weeks, months and years to come, I can contribute some useful Midwest Scenery and Airport Data.

=== Airport Layouts ===
So far, I have updated and submitted airport layouts to the [[FlightGear]] project for
*KEQA
*3AU
*KRCP
I have a few more in the works and plan to eventually add terminals and structures to as many Kansas Airports as I have the motivation to finish.

[[File:KEQA_Airport.jpg|350px|KEQA]]

=== Kansas Custom Scenery ===
I also have custom scenery generated from a combination of
*NLCD (land cover)
*OSM (roads, rivers and rail)
*SRTM (water data)
*New CS Data
*custom textures
for the entire state of [https://1drv.ms/u/s!AmL-8HrJKQ8ciRZl-sKR0HWnKPJQ?e=YuJoQ4 Kansas]

This is currently very heavy data but works surprisingly well on my mid level CPU and graphics card.

<gallery caption="Kansas Scenery" mode="packed-hover">
Kansas Scenery One.jpg|ElDorado, Kansas
Kansas Scenery Four.jpg|ElDorado Lake
Kansas Scenery Two.jpg|Kansas City
KICT.jpg|KICT - Wichita Dwight D. Eisenhower National Airport (Mid-Continent)
KIAB.jpg|KIAB - McConnell AFB, Wichita
Kansas.jpg|Kansas Sunset
KMCI.jpg|KMCI - Kansas City International Airport
ElDorado Lake Sunset.jpg|Sunset at ElDorado Lake
Cub Over Eldorado Lake.jpg|Cub at ElDorado Lake
</gallery>

=== Ohio Custom Scenery ===
Custom scenery generated from a combination of
*CS Data based on NLCD (land cover)
*OSM (roads, rivers and rail)
*SRTM (water data to come)
*custom textures
for the entire state of [https://1drv.ms/u/s!AmL-8HrJKQ8ciAqGzBvP5O-kXZRz?e=1i4xWm Ohio].

<gallery caption="Ohio Scenery" mode="packed-hover">
Cincinnati.jpg|KLUK - Cincinnati Municipal Airport
Toledo-Harbor View.jpg|Toledo - Harbor View
Lake_Erie_2.jpg|KCLE Cleveland-Hopkins International Airport
KPCW.jpg|KPCW - Carl R Keller Field
South Lorain.jpg|KCAK Akron-Canton Regional Airport
KBKL.jpg|KSKY - Griffing Sandusky Airport
Lake Erie.jpg|KCMH Port Columbus International Airport
KDAY_James_M._Cox_Dayton_International_Airport.jpg|KDAY James M. Cox Dayton International Airport
KYNG Youngstown-Warren Regional Airport.jpg|KYNG Youngstown - Warren Regional Airport
Near Youngstown.jpg|Near Youngstown
</gallery>

== Addon Development ==
[[File:Cargo-position.jpg|thumb|501x501px|Cargo placement UI.]]

=== Cargo Towing Addon ===
Cargo Towing is now in the form of an addon. The Cargo Towing addon adds the ability for any helicopter to haul cargo from either a slung position. The AirCrane retains a docked cargo capability. Now included, hooks for adding load weight and winched rope length to any aircraft that has the instrumentation do display this data.. Using the Wildfire addon in addition to Cargo Towing addon, you can use a bambi bucket for firefighting. Also added a special feature that allows the AirCrane to load water from any of the 3 Water Tank cargos. The addon extended the sling load Cargo Towing capability to any aircraft.

[https://svn.code.sf.net/p/flightgear/fgaddon/trunk/Addons/CargoTowingAddon Cargo Towing Addon]

[[File:Intercept1.jpg|left|thumb|357x357px]]

=== Mission Generator Addon ===
The Mission Generator addon is a continuing project that is designed to give FlightGear some ability to easily generate different mission orientated challenges. So far there are currently ...

* Wildfire events. Using Wildfire and Cargo Towing addons you have the capability to put these fire out using various aircraft and methods. Including Water drop with tankers, bambi bucket or the Aircrane's jet spray boom apparatus.
* Search and rescue, water and land, including life flight to hospital.
* Intercept missions, including an airspace intrusion

All missions are fully customizable and can be either on demand or set a random range and time for the event, per event type.

[https://svn.code.sf.net/p/flightgear/fgaddon/trunk/Addons/MissionGenerator Mission Generator Addon]

[[File:Intercept3.jpg|thumb|400x400px|Vectored to target.|left]]
[[File:Intercept2.jpg|left|thumb|408x408px|Target data vectoring.]]

[[File:Misson Generator Help.png|thumb|Help visual for the Mission Generator UI|border|center|820x820px]]
=== WildfireAddon ===
[[File:Bambi-bucket-full.jpg|thumb|761x761px|Bucket fills by submerging in any water source.]]The Wildfire addon originally was a FlightGear stand alone nasal module written and developed by Anders Gidenstam and incorporated in fgdata

Credit:

Wildfire - A cellular automaton forest fire model with the ability to spread over the multiplayer network.Written and developed by Anders Gidenstam (anders(at)gidenstam.org) Copyright (C) 2007 - 2012 Anders Gidenstam (anders(at)gidenstam.org) Modifications and conversion to Addon by Wayne Bragg wlbragg<blockquote>

The cellular automata model used here is loosely based on A. Hernandez Encinas, L. Hernandez Encinas, S. Hoya White, A. Martin del Rey, G. Rodriguez Sanchez, "Simulation of forest fire fronts using cellular automata", Advances in Engineering Software 38 (2007), pp. 372-378, Elsevier.</blockquote>[https://svn.code.sf.net/p/flightgear/fgaddon/trunk/Addons/Wildfire Wildfire Addon]
[[File:Bambi-bucket-dump.jpg|thumb|455x455px|Dumping the bucket will suppress the wildfire..|alt=|none]]

[[File:Bambi-busket2.jpg|thumb|741x741px|View from engineering station.|alt=|none]]

== Future Endeavors ==

=== Surprise ===

[[File:Iss4.jpg|512px|International Space Stay|alt=|left|thumb]][[File:Aircrane Over the Mojave.jpg|1024px|AirCrane over the Mojave Desert in the Palmdale/Landcaster area. Compositor shadows, OSM2CITY models, ALS lighting and Orthographics texture.]]

[[File:Aircrane Over the Mojave 2.jpg|1024px|AirCrane over the Mojave Desert in the Palmdale/Landcaster area. Compositor shadows, OSM2CITY models, ALS lighting and Orthographics texture.]]
{{-}}

Watch this page for future releases of various project for [[FlightGear]]. I hope to have a couple of surprises just around the bend.

2024-10-23T17:40:37Z

Wlbragg:

{{Note|The [[Canvas]] system is a more recent, and much simpler, option to add immatriculation support to an aircraft using CanvasText/osgText nodes (typically, 5 lines of code), for details see: [[Howto:Dynamic Liveries via Canvas]]}}

[[File:DR400_immat.jpg|thumb|270px|Immatriclation at the tail of the [[Robin DR400|DR400]].]]

'''Immatriculation''' makes it possible to add a user changeable callsign to an [[aircraft]].

Examples are based on the [[Robin DR400]].

<font color="red">'''Please note that this article does not cover the model part! Following only this page, will not result in a working immatriclation!'''</font>

=== -set.xml ===
Add this to enable the immatriculation:
<model>
<immat>true</immat>
</model>

Let's make a menu button:
<menubar>
<default>
<menu n="100">
<label>Robin DR400</label>
<enabled type="bool">true</enabled>
<item>
<label>Immatriculation</label>
<binding>
<command>nasal</command>
<script>dr400.immat_dialog.toggle()</script>
</binding>
</item>
</menu>
</default>
</menubar>

</sim>

And place a path to the immat.nas file we will create in the next step.
<nasal>
<file>Aircraft/DR400/Nasal/immat.nas</file>
</nasal>

=== immat.nas ===
<pre>
# ===========================
# Immatriculation by Zakharov
# ===========================

var refresh_immat = func {
var immat = props.globals.getNode("/sim/model/immat",1).getValue();
var immat_size = size(immat);
if (immat_size != 0) immat = string.uc(immat);
for (var i = 0; i < 6; i += 1) {
if (i >= immat_size)
glyph = -1;
elsif (string.isupper(immat[i]))
glyph = immat[i] - `A`;
elsif (string.isdigit(immat[i]))
glyph = immat[i] - `0` + 26;
else
glyph = 36;
props.globals.getNode("/sim/multiplay/generic/int["~i~"]", 1).setValue(glyph+1);
}
}

var immat_dialog = gui.Dialog.new("/sim/gui/dialogs/DR400/status/dialog",
"Aircraft/DR400/Dialogs/immat.xml");

setlistener("/sim/signals/fdm-initialized", func {
if (props.globals.getNode("/sim/model/immat") == nil) {
var immat = props.globals.getNode("/sim/model/immat",1);
var callsign = props.globals.getNode("/sim/multiplay/callsign").getValue();
if (callsign != "callsign") immat.setValue(callsign);
else immat.setValue("F-GHYQ");
}
refresh_immat();
setlistener("sim/model/immat", refresh_immat, 0);
},0);
</pre>

<font color="red">'''Example of the model part mentioned above!'''</font>
There are four areas of the aircraft registration that are covered in this file, the Wing, left and right fuselage, and interior panel.
The list of model parts on the left side of the image are the individual letter or number positions of the final registration identification string.
[[File:Immatriculation_modeling_in_Blender.png|thumb|640px|Immatriculation Blender Model For the DR400.]]

[[Category:Aircraft enhancement]]
[[Category:Nasal]]

Immatriculation

2024-10-23T17:39:41Z

Wlbragg:

{{Note|The [[Canvas]] system is a more recent, and much simpler, option to add immatriculation support to an aircraft using CanvasText/osgText nodes (typically, 5 lines of code), for details see: [[Howto:Dynamic Liveries via Canvas]]}}

[[File:DR400_immat.jpg|thumb|270px|Immatriclation at the tail of the [[Robin DR400|DR400]].]]

'''Immatriculation''' makes it possible to add a user changeable callsign to an [[aircraft]].

Examples are based on the [[Robin DR400]].

<font color="red">'''Please note that this article does not cover the model part! Following only this page, will not result in a working immatriclation!'''</font>

=== -set.xml ===
Add this to enable the immatriculation:
<model>
<immat>true</immat>
</model>

Let's make a menu button:
<menubar>
<default>
<menu n="100">
<label>Robin DR400</label>
<enabled type="bool">true</enabled>
<item>
<label>Immatriculation</label>
<binding>
<command>nasal</command>
<script>dr400.immat_dialog.toggle()</script>
</binding>
</item>
</menu>
</default>
</menubar>

</sim>

And place a path to the immat.nas file we will create in the next step.
<nasal>
<file>Aircraft/DR400/Nasal/immat.nas</file>
</nasal>

=== immat.nas ===
<pre>
# ===========================
# Immatriculation by Zakharov
# ===========================

var refresh_immat = func {
var immat = props.globals.getNode("/sim/model/immat",1).getValue();
var immat_size = size(immat);
if (immat_size != 0) immat = string.uc(immat);
for (var i = 0; i < 6; i += 1) {
if (i >= immat_size)
glyph = -1;
elsif (string.isupper(immat[i]))
glyph = immat[i] - `A`;
elsif (string.isdigit(immat[i]))
glyph = immat[i] - `0` + 26;
else
glyph = 36;
props.globals.getNode("/sim/multiplay/generic/int["~i~"]", 1).setValue(glyph+1);
}
}

var immat_dialog = gui.Dialog.new("/sim/gui/dialogs/DR400/status/dialog",
"Aircraft/DR400/Dialogs/immat.xml");

setlistener("/sim/signals/fdm-initialized", func {
if (props.globals.getNode("/sim/model/immat") == nil) {
var immat = props.globals.getNode("/sim/model/immat",1);
var callsign = props.globals.getNode("/sim/multiplay/callsign").getValue();
if (callsign != "callsign") immat.setValue(callsign);
else immat.setValue("F-GHYQ");
}
refresh_immat();
setlistener("sim/model/immat", refresh_immat, 0);
},0);
</pre>

<font color="red">'''Example of the model part mentioned above!'''</font>
There are four areas of the aircraft registration that are covered in this file, the Wing, left and right fuselage, and interior panel.
The list of model parts on the left side of the image are the individual letter or number positions of the final registration identification string.
[[File:Immatriculation_modeling_in_Blender.png|thumb|640px|Immatriculation Blender Model For the c172p.]]

[[Category:Aircraft enhancement]]
[[Category:Nasal]]

Immatriculation

2024-10-23T17:38:47Z

Wlbragg:

File:Immatriculation modeling in Blender.png

2024-10-23T17:36:17Z

Wlbragg: Wlbragg uploaded a new version of File:Immatriculation modeling in Blender.png

=={{int:filedesc}}==
{{Information
|description={{en|1=Showing the immatriculation modeling in Blender that is used in on the c172p model.}}
|date=2024-10-23
|source={{own}}
|author=[[User:Wlbragg|Wlbragg]]
|permission=
|other versions=
}}

=={{int:license-header}}==
{{self|cc-by-sa-4.0}}

[[Category:Immatriculation modeling Blender registration]]

Immatriculation

2024-10-23T16:21:40Z

Wlbragg:

{{Note|The [[Canvas]] system is a more recent, and much simpler, option to add immatriculation support to an aircraft using CanvasText/osgText nodes (typically, 5 lines of code), for details see: [[Howto:Dynamic Liveries via Canvas]]}}

[[File:DR400_immat.jpg|thumb|270px|Immatriclation at the tail of the [[Robin DR400|DR400]].]]

'''Immatriculation''' makes it possible to add a user changeable callsign to an [[aircraft]].

Examples are based on the [[Robin DR400]].

<font color="red">'''Please note that this article does not cover the model part! Following only this page, will not result in a working immatriclation!'''</font>

=== -set.xml ===
Add this to enable the immatriculation:
<model>
<immat>true</immat>
</model>

Let's make a menu button:
<menubar>
<default>
<menu n="100">
<label>Robin DR400</label>
<enabled type="bool">true</enabled>
<item>
<label>Immatriculation</label>
<binding>
<command>nasal</command>
<script>dr400.immat_dialog.toggle()</script>
</binding>
</item>
</menu>
</default>
</menubar>

</sim>

And place a path to the immat.nas file we will create in the next step.
<nasal>
<file>Aircraft/DR400/Nasal/immat.nas</file>
</nasal>

=== immat.nas ===
<pre>
# ===========================
# Immatriculation by Zakharov
# ===========================

var refresh_immat = func {
var immat = props.globals.getNode("/sim/model/immat",1).getValue();
var immat_size = size(immat);
if (immat_size != 0) immat = string.uc(immat);
for (var i = 0; i < 6; i += 1) {
if (i >= immat_size)
glyph = -1;
elsif (string.isupper(immat[i]))
glyph = immat[i] - `A`;
elsif (string.isdigit(immat[i]))
glyph = immat[i] - `0` + 26;
else
glyph = 36;
props.globals.getNode("/sim/multiplay/generic/int["~i~"]", 1).setValue(glyph+1);
}
}

var immat_dialog = gui.Dialog.new("/sim/gui/dialogs/DR400/status/dialog",
"Aircraft/DR400/Dialogs/immat.xml");

setlistener("/sim/signals/fdm-initialized", func {
if (props.globals.getNode("/sim/model/immat") == nil) {
var immat = props.globals.getNode("/sim/model/immat",1);
var callsign = props.globals.getNode("/sim/multiplay/callsign").getValue();
if (callsign != "callsign") immat.setValue(callsign);
else immat.setValue("F-GHYQ");
}
refresh_immat();
setlistener("sim/model/immat", refresh_immat, 0);
},0);
</pre>

<font color="red">'''Example of the model part mentioned above!'''</font>
There are four areas of the aircraft registration that are covered in this file, the Wing, left and right fuselage, and interior panel.
The list of model parts on the left side of the image are the individual letter or number positions of the final registration identification string.
[[File:Immatriculation_modeling_in_Blender.png|thumb|640px|Immatriculation Blender Model For the c172P.]]

[[Category:Aircraft enhancement]]
[[Category:Nasal]]

File:Immatriculation modeling in Blender.png

2024-10-23T15:54:52Z

Wlbragg: Uploaded own work with UploadWizard

Howto:Create WS3.0 terrain

2024-08-01T03:16:03Z

Wlbragg: /* Optional Python script to process NLCD for the USA */

{{WS30 Navbar}}
This article provides instructions on how to generate basic WS3.0 terrain.

WS3.0 terrain consists of three parts:

# A terrain mesh consisting of a landclass texture draped over an elevation model.
# A high resolution water raster used to show water features such as rivers, lakes and coastline with more definition
# Line features such as roads and railways.

The terrain is generated by a set of tools that are packaged in a docker image for convenience.[[File:Diagram-export-21-12-2023-16 29 37.png|thumb|Basic WS3.0 Scenery Generation Process]]

== Prerequisites ==

=== Set up a Workspace ===
Create a directory with the following sub-directories:

* <code>cache</code>
* <code>data</code>
* <code>output/vpb</code>
* <code>output/Terrain</code>

=== Docker ===

# Install [https://docs.docker.com/get-started/ Docker] on your platform.
#Pull the docker image by running the following command

docker pull flightgear/ws30-vpb-generator:latest
Optionally, if you are hitting rate limits:
#Create an account on https://hub.docker.com/. (Note that you will need to click on an email verification link before you can log in for the first time)
#Run <code>docker login</code> before the '''docker pull''' command above

== Getting the base data ==
You need two pieces of data for the area of scenery you are generating:

# An elevation model (aka DEM). This indicates what altitude each point of the surface is.
# Landclass data showing what type of terrain is at each point of the surface. This is often either a Raster (effectively a texture), or vector data.

=== Elevation Model ===
Download the NASADEM elevation model for the area of scenery you wish to generate. This is available in 1x1 degree blocks from [https://lpdaac.usgs.gov/products/nasadem_hgtv001/ here], and with an interactive browser [https://search.earthdata.nasa.gov/search here].

Unzip the files into the <code>data</code> directory.

=== Landclass Raster ===
Download an landclass raster for the area of scenery you wish to generate.

* For Europe, use of [https://land.copernicus.eu/pan-european/corine-land-cover/clc2018 CORINE] is recommended.
* For the USA [https://www.mrlc.gov/viewer/ NLCD] is recommended
* Sentinel-2 data is available for the entire world via [https://livingatlas.arcgis.com/landcoverexplorer/ ESRI], but has limited set of landclasses.

Put these into the <code>data</code> directory.

More detailed terrain can be created by modifying the landclass raster, and/or generating a new raster from vector data. These processes are discussed below.

== Generating Terrain ==
To generate terrain you need to run the tools within the docker container we installed above. The docker image is like a small, independent virtual computing environment running within your system. This particular docker image has all the scenery generation tools already installed.

=== Running the docker container ===
Firstly, get the container running from the directory containing your <code>cache</code>, <code>data</code> and <code>output</code> directories:
docker run --rm --mount "type=bind,source=`pwd`/data,target=/home/flightgear/data,readonly" --mount "type=bind,source=`pwd`/output,target=/home/flightgear/output" --mount "type=bind,source=`pwd`/cache,target=/home/flightgear/cache" -it flightgear/ws30-vpb-generator:latest /bin/bash

You should now find yourself in a bash shell within your container. You should see data and output directories which are linked to the directories you created earlier:
flightgear@ddcac77f7d5e:~$ ls
cache data output bin scripts

=== Building the terrain ===
To build the terrain mesh, use the <code>genVPB.py</code> tool from inside the docker container.
Usage: genVPB.py <min-lon> <min-lat> <max-lon> <max-lat> <input-raster> [--reclass <reclass>] [--coastline <coastline>] [--hgt-dir <hgt-dir>] [--output-dir <output-dir>]
<min-lon> <min-lat> <max-lon> <max-lat> - bounding box of scenery to generated
<input-raster> - Input landclass raster
[coastline] - Optional coastline polygon data (.osm) to clip against. E.g. from <nowiki>https://osmdata.openstreetmap.de/data/land-polygons.html</nowiki>
[reclass] - Optional file containing a set of rules for reclassifying the raster similar to r.reclass. See fgmeta/ws30/mappings/
[hgt-dir] - Optional directory containing unzipped NASADEM HGT files. Defaults to '/home/flightgear/data'.
[output-dir] - Optional directory for output scenery. Defaults to '/home/flightgear/output/vpb'

For example, to generate a piece of terrain around Edinburgh (latitude 55.5, longitude 3 degrees West)
./scripts/genVPB.py -3 55 -2 56 ./data/uk_wgs84_10m_N54.tif
If you are using anything other than a CORINE raster you will need to reclassify the data to match the landclasses used by FlightGear. Those classes are defined in [https://sourceforge.net/p/flightgear/fgdata/ci/next/tree/Materials/base/landclass-mapping.xml Materials/base/landclass-mapping.xml]. You can reclassify them using the files in the [https://sourceforge.net/p/flightgear/fgmeta/ci/next/tree/ws30/mappings/ scripts/mappings] directory. E.g. to reclassify NLCD2019 data you can use <code>--reclassify ./scripts/mappings/nlcd2019.txt</code>.

genVPB.py will output the data to the output/vpb directory, in which you should find a series of files and directories.

=== Adding water ===
The terrain mesh does not have highly detailed water features - as typically the source data has a resolution of 10-25m. Water features are generated from OpenStreetMap data. To generate water features simply run the genwaterraster.py command.
Usage: ./scripts/genwaterraster.py <scenery_dir> <lon1> <lat1> <lon2> <lat2> [Cache]
<scenery_dir> Scenery directory to write to
<lon1> <lat1> Bottom left lon/lat of bounding box
<lon2> <lat2> Top right lon/lat of bounding box
[Cache] Optional directory to cache OpenStreetMap data to reduce API requests
The scenery directory should be your output/vpb directory created above. For example, to generate water for the area around Edinburgh:
./scripts/genwaterraster.py ./output/vpb -4 55 -3 56 ./cache
Inside the ./output/vpb directory there should be a set of directories and .png files.

=== Adding roads and railways ===
The terrain mesh does not have any line features - things like roads. These are generated separately from OpenStreetMap data. To generate line features simply run the genroads.py command:
Usage: ./scripts/genroads.py <scenery_dir> <lon1> <lat1> <lon2> <lat2> [Cache]
<scenery_dir> Scenery directory to write to
<lon1> <lat1> Bottom left lon/lat of bounding box
<lon2> <lat2> Top right lon/lat of bounding box
[Cache] Optional directory to cache OpenStreetMap data to reduce API requests

The scenery directory should be your output/Terrain directory created above. For example, to generate roads for the area around Edinburgh:
./scripts/genroads.py ./output/Terrain -4 55 -3 56 ./cache
Inside the ./output/Terrain directory there should be a set of directories and, .STG files text files.

==Running FlightGear with the new WS3.0 Terrain==
To test the new terrain, simply include the output directory in your scenery path and run FlightGear with the <code>--prop:/scenery/use-vpb=true</code> to enable WS3.0.

== Advanced Techniques ==
The following sections describe more complex techniques to generate higher quality WS3.0 terrain. Almost all of them involve using different data sources to generate a more detailed landclass raster before running the final scenery generation processes described above. Generating a highly detailed landclass raster is where the magic happens.

Most techniques use gdal or grass to modify the raster/vector data, typically using the QGIS program.

=== Using a different elevation model ===
If you are using another elevation model other than NASAEM, then you may need to re-project it using QGIS/gdalwarp to the WGS84 CRS (aka EPSG:4326).

=== Landclass Data Requirements ===
For any landclass data we need to ensure the data is in the correct format. That means:

# Is a Raster (geotiff) rather than Vector data. This raster will become the texture on the terrain that the terrain shaders do their magic on.
# Uses the WGS84 Coordinate Reference System. The ensures that the terrain generation step is efficient.
# Has the correct landclass values for each terrain type. We use a set of values based on the CORINE raster set, defined in [https://sourceforge.net/p/flightgear/fgdata/ci/next/tree/Materials/base/landclass-mapping.xml Materials/base/landclass-mapping.xml].

Below is a quick table showing what steps you need to take for common landclass data sources.
{| class="wikitable"
|+
!Landclass Data
!Warp to WGS84 required?
!Landclass re-classification Required?
!Raster Simplification Required?
!Conversion to Raster Required?
|-
|CORINE Raster
|Yes
|No
|No
|No
|-
|CORINE Vector
|Yes
|Yes
|No
|Yes
|-
|NLCD
|No
|Yes
|Yes
|No
|-
|Sentinel-2
|Yes
|Yes
|No
|No
|}
Conversion to Raster must be done manually. Converting to WGS84 and the correct landclasses ''can'' be done by the genVPB.py script, but slows down scenery generation. Therefore if you are planning to generate scenery multiple times it is best to pre-process the files yourself.

The easiest way to do these operations is using QGIS, which is available for most platforms. If you are scripting a toolchain, the QGIS tools include command-line equivalents for all commands.

When using QGIS, set the Project CRS to WGS84 (aka EPSG:4326). You can then add layers of Raster or Vector data from files from the <code>Layer->Add Layer</code> menu. When performing any operations, <u>always</u> write out the data to a real file so you can go back to it later. Disk space is cheap :).

=== Warping Raster Layers to WGS84 ===
genVPB.py will automatically convert to WGS84 if it detects a landclass raster with a different Coordinate System. However, you may wish to do so manually instead so all your data is on the same coordinate system.

In QGIS you can warp a raster layer to a different CRS using the Raster->Projections-Warp (Reproject) tool.

Select the following options in the dialog:

* Input Layer - Check you have selected the correct layer. The CRS is shown at the right.
* Target CRS - set to <code>EPSG:4326 - WGS84</code>, which should be your project CRS.
* Resampling Method to Use - Nearest Neighbour. (Landclass data is not like normal images. You don't want to interpolate between values.)
* Nodata value for output bands - 0.0 (This means that any data at the edges will be Ocean, usually a reasonable default)
* Advanced Parameters - No Compression
* Output data type - Byte
* Reprojected - Choose a target filename

Alternatively you can do this step from the commandline.
gdalwarp -t_srs EPSG:4326 -dstnodata 0.0 -r near -ot Byte -of GTiff -co COMPRESS=NONE -co BIGTIFF=IF_NEEDED /home/stuart/FlightGear/VPB/data/CORINE/u2018_clc2018_v2020_20u1_raster100m/DATA/U2018_CLC2018_V2020_20u1.tif /home/stuart/FlightGear/VPB/data/scratch/corine_WGS84.tif

=== Warping Vector Layers to WGS84 ===
genVPB.py will automatically convert to WGS84 if it detects a landclass raster with a different Coordinate System. However, you may wish to do so manually instead so all your data is on the same coordinate system.

In QGIS you can warp a vector layer using Vector->Data Management Tools->Reproject Layer.

Set the following options in the dialog:

* Input Layer - Check you have selected the correct layer. The CRS is shown at the right.
* Target CRS - set to <code>EPSG:4326 - WGS84</code>, which should be your project CRS.
* Reprojected - Choose a target filename

=== Reclassifying Vector Layers ===
For CORINE vector data in particular, the attributes used in the vector data are not the same as those used by the CORINE Raster data. So we need to create a new attribute on the data.
[[File:Field Calculator.png|thumb|QGIS Field Calculator]]
To do this

* Select <code>Layer->Open Attribute Table</code>. You should see a table with multiple columns. Each row represents a feature in the data.
* Click on <code>Toggle Editing Mode</code> (Ctrl + E), on the top left of the dialog
* Click on <code>New Field</code> (Ctrl+W) and create a new field called Landclass of type "Whole Number (Integer)". This will create a new column which we will populate with the correct landclass data
* Click on the <code>Open Field Calculator</code> button (Ctrl + I). (If you get an error about only being able to create Virtual fields, go back to the Layer menu, export it and open the exported file).
* Select the following options:
** Update Existing Field
** Select the Landclass field you just created.
** Copy the contents of https://sourceforge.net/p/flightgear/fgmeta/ci/next/tree/ws30/mappings/corine_vector.txt into the Expression box (without the comment lines starting with <code>#</code>). This is just some simple code to set the attribute correctly. The code should be correct for CORINE vector data. If your data is from other sources you will need to work out how you want to map your source data landclasses to the CORINE ones. [https://sourceforge.net/p/flightgear/fgdata/ci/next/tree/Materials/base/landclass-mapping.xml Materials/base/landclass-mapping.xml] can be used as a guide.
* Select <code>OK</code>. You should see that your landclass column is now populated with the landclass data.
* Select <code>Layer->Save Layer Edits</code> to save you changes

=== Creating a Raster from a Vector Layer ===
To create a Raster from a Vector Layer select <code>Raster->Conversion->Rasterize (Vector to Raster)</code>.
[[File:QGIS Rasterize (Vector to Raster).png|thumb|Creating a Raster from a Vector Layer - QGIS Rasterize]]
Select the following options:

* Input Layer - correct layer, check CRS
* Field to use for burn-in value - select the Landclass column you created above.
* Output raster size units. This is going to set the resolution of your raster. You can work out the resolution in two different ways:
** Select "Georeferenced units" and determine how many degrees each pixel is in latitude and longitude.
** Select "Pixels" and determine the size of raster you want in pixels. [https://www.nhc.noaa.gov/gccalc.shtml This] is a good calculator to help.
* Width/Horizontal Resolution. Enter the values you've calculated for the horizontal resolution (longitudinal), or the width of the raster
* Height/Vertical Resolution. Enter the values you've calculated for the vertical resolution (latitude or the height of the raster
* Output extent - Select an option from the box on the right. You can edit the text afterwards. Best practise is to create long thin strips of 1 degree latitude in height, as this makes subsequent processing much easier.
* Assign a specific nodata value to output bands - Select 0.0 for Ocean. CORINE vector data in particular has a lot of nodata for Oceans
* Advanced Parameters - No Compression
* Output data type - Byte
* Rasterized - Select a new filename

=== Simplifying a Raster Layer ===
Some Raster Landclass data (NLCD included) has too much noise - in particular large US highway systems are identified as Urban areas.

To smooth it out we can use the GRASS <code>n.neighbors</code> function from the Processing Toolbox in QGIS.

Select the following options:

* Input Layer - correct layer, check CRS
* Neighborhood operation - median. (This is not a normal image, so using an average will result in weird values)
* Neighborhood size - 5.
* Neighbors - Select a new filename.

=== Clipping a Raster Layer with OSM Data for Land (Corine) ===
The Corine dataset does not match OSM coastlines exactly. The following multi-stage process makes sure, that no Corine land-use is in the water as defined by OSM.

==== Download OSM Land Data ====

Download land polygons based on OSM data as a Shapefile from [https://osmdata.openstreetmap.de/data/land-polygons.html Land Polygons] and make sure to pick the WGS84 projected download with split polygons ("Large polygons are split, use for larger scales"). Once downloaded unzip the content into a directory.

==== Reclassifying the OSM Land Data Vector Layer ====
To do this:

* Select <code>Layer->Open Attribute Table</code>. You should see a table with multiple columns. Each row represents a feature in the data.
* Click on <code>Toggle Editing Mode</code> (Ctrl + E), on the top left of the dialogue
* Click on <code>New Field</code> (Ctrl+W) and create a new field called Landclass of type <code>Integer (32 bit)</code>. This will create a new column which we will populate with the correct land class data
* On top of the table on the left side choose "Landclass" in the drop-down menu, then input <code>2</code> into the field to the right and then press button "Update" all to the left of this field.
* Wait a bit and the close the dialogue.
* Select <code>Layer->Save Layer Edits</code> to save you changes

==== Convert the Land Data from Vector to Raster ====
Do the same as in chapter "Creating a Raster from a Vector Layer" above. The only difference is that the Input layer will be the land data polygons and you need to choose a different file name for the "Rasterized" (e.g. osm_land_scotland_5m.tif)

==== Remove Novalue Entries in the Land Data Raster ====
To do this:

* Select <code>Processing->Toolbox</code>. You should see a new box on the right side.
* Write "gdal_calc" in the search box and you should see an entry "Raster calculator". Double click on it and you will get a new dialogue window.
* In this dialogue:
** For "Input layer A" choose the raster from the previous chapter ((e.g. osm_land_scotland_5m.tif)
** In field "Calculation in gdalnumeric ..." write: <code>greater(A,0) * A</code>
** In field "Output raster type" choose <code>Byte</code>
** In field "Advanced Parameters" choose Profile <code>No compression</code>
** In field "Additional command-line parameters" write: <code>--hideNoData</code>
** In field "Calculated" choose a file (e.g. osm_land_scotland_allvalues_5m.tif)
** (In the "GDAL/OGR console call" it will have something similar to the follwing - just with different paths: <code>gdal_calc.py --overwrite --calc "greater(A ,0) * A" --format GTiff --type Byte -A /home/vanosten/custom-fg-scenery/data/osm_land_scotland_5m.tif --A_band 1 --co COMPRESS=NONE --co BIGTIFF=IF_NEEDED --hideNoData --outfile /home/vanosten/custom-fg-scenery/data/osm_land_scotland_allvalues_5m.tif</code>
** Press the "Run" button - and when complete close the dialogue.

You should now see a map only black and white. You can check for correctness by pressing <code>CTRL+SHIFT+I</code> to get a cursor with an arrow and an "i". First make sure the new raster is selected on the left side. Next click on the sea/ocean and then check in the "Identify Results" window on the right that the value is < 2. The click on the land and check that the value is 2.

==== Create the Final Clipped Corine Raster Against OSM Land Data =====
Do the following:

* In QGIS make sure that you have only the following two layers: the basis Corine raster (see chapter "Creating a Raster from a Vector Layer" - here e.g. corine_raster_scotland_5m.tif) and plus the raster from the previous step (e.g. osm_land_scotland_allvalues_5m.tif)
* Select <code>Raster->Raster Calculator ...</code> and a corresponding dialogue will open showing on the left hand side the two rasters.
* Choose a new "Output layer" (e.g. corine_raster_scotland_clipped_5m.tif).
* In the "Raster Calculator Expression" field input: <code>if ("osm_land_scotland_all_data_5m@1" < 2, 44, "corine_raster_scotland_5m@1")</code>
* Press button "OK" and wait a while (you will see a new dialogue with showing the progress.

Done. You have now a raster (e.g. corine_raster_scotland_clipped_5m.tif) which does not have land in areas, where OSM data has sea/ocean.

=== Reclassifying a Raster Layer ===
WS3.0 uses CORINE landclass values. If using data from other sources it needs to be reclassified to the correct values. genVPB.py has an option to do this, but you may wish to do so manually.

To do this select <code>GRASS->Raster->r.reclass</code> from the Processing Toolbox.

Select the following options:

* Input Raster Layer - correct layer, check CRS
* Reclass rules text - copy in the contents of https://sourceforge.net/p/flightgear/fgmeta/ci/next/tree/ws30/mappings/nlcd2019.txt. Or an appropriate mapping from your landclass data to CORINE. Note that you can also reference a file using the "File containing reclass rules" option. Note a mapping of 22 24 = 1 is the same as 22 and 24 = 1. For a range of 22 to 24 use 22 23 24 = 1.
* Reclassified - Select a new filename.

(If this doesn't work a similar function is available in the Processing Toolbox under <code>Raster analysis->Reclassify by table</code>. However this doesn't save your table once you close the dialog, and entries have to be manually entered individually which takes a lot of effort)

=== Step By Step Procedure for Processing NLCD for the USA using the Raster Calculator, Up sampling and GRASS r.neighbors===

We will use a predetermined file naming convention throughout this procedure for simplicity. You can use a naming convention of your choice.

Obtain your NLCD_2019_Land_Cover from source of choice, https://www.mrlc.gov/viewer/

NLCD_2019_Land_Cover data = NLCD_2019_Land_Cover_California-Southern.tiff

Obtain your NLCD_2016_Tree_Canopy from source of choice, https://www.mrlc.gov/viewer/

NLCD_2016_Tree_Canopy data = NLCD_2016_Tree_Canopy_California-Southern.tiff

Make Deciduous Layer 255 to 0 EPSG:5070 - NAD83 - Raster Calculator

("NLCD_2016_Tree_Canopy_California-Southern@1" = 255) * 0 + ("NLCD_2016_Tree_Canopy_California-Southern@1" != 255) * "NLCD_2016_Tree_Canopy_California-Southern@1"

Output CRS = EPSG:5070 - NAD83

Output Layer = California-Southern_deciduous-coast-clipped.tiff

OK

Make Tree Canopy Layer EPSG:5070 - NAD83 - Raster Calculator

Use the prevailing tree type, for Deciduous use * 41, for Evergreen use * 42, for MixedForest use 43.

Later we add this Tree Canopy layer to the raster where there is currently no trees in the NLCD

If ocean bordering with index 255 -

("California-Southern_deciduous-coast-clipped@1" != 0) * 42

Else

("NLCD_2016_Tree_Canopy_California-Southern@1" != 0) * 42

Output CRS = EPSG:5070 - NAD83

Output Layer = California-Southern_deciduous.tiff

OK

Make Clipped Ocean Frontage EPSG:5070 - NAD83 - GIMP

Source = NLCD_2019_Land_Cover_California-Southern.tiff

Flood fill index 11 (blue) with 0 (black)

Export tiff = California-Southern_coast-clipped.tiff

Warp NLCD Deciduous and Base NLCD from EPSG:5070 - NAD83 / Conus Albers NLCD to Project CRS: EPSG:4326 - WGS 84 - WARP

If ocean bordering with index 255

Warp California-Southern_coast-clipped.tiff

Else

Warp NLCD_2016_Tree_Canopy_California-Southern.tiff

and

Warp Deciduous NLCD from EPSG:5070 - NAD83 / Conus Albers NLCD to Project CRS: EPSG:4326 - WGS 84

Warp California-Southern_deciduous.tiff

Raster -> Projection - > Warp

Target CRS = EPSG:4326 - WGS 84

Resampling method to use = Nearest Neighbour

Output data type = byte

Reprojected = California-Southern_4326-84.tiff and California-Southern_deciduous_4326-84.tiff

Run

Combine "California-Southern_deciduous_4326-84.@1" and "California-Southern_4326-84@1" - Raster Calculator

This step adds the Tree Canopy layer to the existing NLCD raster only where there is currently no existing tree land cover.

(("California-Southern_4326-84@1" > 0 AND "California-Southern_4326-84@1" != 41 AND "California-Southern_4326-84@1" != 42 AND "California-Southern_4326-84@1" != 43) AND "California-Southern_deciduous_4326-84@1" > 0) * "California-Southern_deciduous_4326-84@1" + ("California-Southern_4326-84@1" = 41 OR "California-Southern_4326-84@1" = 42 OR "California-Southern_4326-84@1" = 43 OR "California-Southern_deciduous_4326-84@1" <= 0) * "California-Southern_4326-84@1"

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_adjusted.tiff

OK

Re-class urban 21, 22, 23, 24 = grass 26 - Raster Calculator

This step is being used to erase all the man made clutter on the raster with grassland. Instead of grassland you can replace it with any type of prevailing groundcover in the area your processing.

For example this places a road easement of grass on all the road clutter that was originally in the raster. If your processing a desert area you might want to use dirt or sand instead of grassland. For Arctic areas maybe use tundra.

Part of the clutter gets replaced later in the process.

("California-Southern_adjusted@1" = 11) * 41 + ("California-Southern_adjusted@1" = 12) * 34 + ("California-Southern_adjusted@1" = 21) * 26 + ("California-Southern_adjusted@1" = 22) * 26 + ("California-Southern_adjusted@1" = 23) * 26 + ("California-Southern_adjusted@1" = 24) * 26 + ("California-Southern_adjusted@1" = 31) * 27 + ("California-Southern_adjusted@1" = 41) * 26 + ("California-Southern_adjusted@1" = 42) * 24 + ("California-Southern_adjusted@1" = 43) * 25 + ("California-Southern_adjusted@1" = 51) * 30 + ("California-Southern_adjusted@1" = 52) * 29 + ("California-Southern_adjusted@1" = 71) * 26 + ("California-Southern_adjusted@1" = 72) * 32 + ("California-Southern_adjusted@1" = 73) * 31 + ("California-Southern_adjusted@1" = 74) * 31 + ("California-Southern_adjusted@1" = 75) * 32 + ("California-Southern_adjusted@1" = 81) * 18 + ("California-Southern_adjusted@1" = 82) * 19 + ("California-Southern_adjusted@1" = 90) * 25 + ("California-Southern_adjusted@1" = 95) * 35

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_reclass-grass.tiff

OK

Make Urban layer

("California-Southern_adjusted@1" = 21) * 10 + ("California-Southern_adjusted@1" = 22) * 1 + ("California-Southern_adjusted@1" = 23) * 1 + ("California-Southern_adjusted@1" = 24) * 2

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_urban.tiff

OK

Remove the road clutter with r.neighbor

GRASS/Raster/r.neighbor

Neighborhood operation = median

Neighborhood size = 7

Neighbors = California-Southern_urban-only.tiff

Run

Combine "California-Southern_reclass-urban@1" and "California-Southern_reclass-grass@1"

This step replaces the "cleaned up" urban areas back into the raster.

("California-Southern_urban-only@1" < 1) * "California-Southern_reclass-grass@1" + ("California-Southern_urban-only@1" != 0) * "California-Southern_urban-only@1"

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_adjusted-combined.tiff

OK

Up sample California-Southern_adjusted-combined.tiff 2x, 4x or 8x resolution 0.00008309125

Active Layer = California-Southern_adjusted-combined.tiff

Layer/Save As

2x Horizontal = 0.00008309125

2x Vertical = 0.00008309125

4x Horizontal = 0.00010637625

4x Vertical = 0.00010637625

8x Horizontal = 0.00005318812

8x Vertical = 0.00005318812

File Name = California-Southern_final-prep-4x.tiff

OK

Simplify California-Southern_4326-84-hd.tiff with r.neighbor

GRASS/Raster/r.neighbor

Neighborhood operation = median

Neighborhood size = 7

Neighbors = California-Southern_urban-only.tiff

Run

Reclass index 0 to 44, leave the rest the same

("California-Southern_4326-84-hd@1" = 0) * 44 + ("California-Southern_4326-84-hd@1" != 0) * "California-Southern_4326-84-hd@1"

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_4326-84-hd-corrected.tiff

OK

<big>Optional HD Fresh Water Option</big>
*Step 1 - Obtain and load hi resolution vector layer
Make sure that the vector layer and the raster layer you will eventually merge to have the same projection. ** Extent can be different if you use the option below.
#Use Top Menu: "Raster" -> "Conversion" -> "Rasterize (vector to raster)" or Processing Toolbox: GDAL -> "Vector conversion" -> "Rasterize (vector to raster)"
# Input layer = Southern-California_water_4326-84
#Fixed value to burn = 41 (water)
# Output raster size units = "Georeferenced units"
#Width/Horizontal resolution = 0.00008309125 (4x the base NLCD resolution )
#Height/Vertical resolution = 0.00008309125 (4x the base NLCD resolution)
#<nowiki>** Output extent = Select the raster layer you will eventually merge with as the "Output extent".</nowiki>
#Under "Rasterized" choose a new filename to save the new raster to. We'll use Southern-California_4326-84-hd-water.tiff.

Click "Run" to generate the hi resolution raster layer.

*Step 2 - Reclass hi res, smoothed water to grass
# Open Raster Calculator
#Under "Output layer" choose a new filename to save the new raster to. We'll use Southern-California_4326-84-hd-nowater.tiff.
#Copy the following into the "Raster Calculator Expression" box. Where Southern-California_4326-84-hd is the name of your hi resolution, smoothed NLCD that includes the water data.
"Southern-California_4326-84-hd@1" * ("Southern-California_4326-84-hd@1" != 41) + 26 * ("Southern-California_4326-84-hd@1 = 41")

Click on "OK". When finished you will have a new raster with the water layer changed to grassland. Another choice would be to change it to sand.

*Step 3 - Combine the hi resolution no water raster and the hi resolution water raster.
# Use Top Menu: "Raster" -> "Miscellaneous" -> "Merge" or Processing Toolbax: GDAL -> Raster Miscellaneous -> Merge
#Input layers = Select "Southern-California_4326-84-hd-nowater@1" and Southern-California_4326-84-hd-water.tiff
#Output data type = "byte".
# Under "Merged" choose a new filename to save the new raster to. We'll use Southern-California_4326-84-hdwater.tiff.

Click "Run" to generate the new merged hi resolution raster layer.

=== Optional Python script to process NLCD for the USA ===
You can optionally process the NLCD by loading this script in the python console of the QGIS Desktop program.

If you save the NLCD tiff's using the same consistent naming convention used in this example it is fairly simple to generate final WS3.0 tiff's to use to generate the scenery.

You will start with an original tree layer and a land cover layer from the MRLC.gov site. See the section "Obtaining NLCD" in the above "Step By Step Procedure for Processing NLCD for the USA using the Raster Calculator" topic.

For example, to start with you will have a couple files named...

NLCD_2021_Tree_Canopy_NewYork.tiff

NLCD_2021_Land_Cover_NewYork.tiff

In the QGIS program the layer names you start with will be...

NLCD_2021_Tree_Canopy_NewYork

NLCD_2021_Land_Cover_NewYork

which will refer to the tiff files.

Running this script will generate a bunch of intermediate files including the last file, that will be the final finished file for running in the VPB script or VPB build environment Docker image.

The last file for this example for New York would be, NLCD_2021_Land_Cover_NewYork_Final-HD_4326.tiff<syntaxhighlight lang="python">
from qgis.core import QgsProject, QgsRasterLayer, QgsCoordinateReferenceSystem, QgsProcessingException
from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry
from osgeo import gdal, osr
import os
import numpy
import numpy as np
import processing
import subprocess

# Define input layer names, change according to your file names
year = '2021';
state = 'Kentucky';
part = '';

tree_canopy_layer_name = 'NLCD_' + year + '_Tree_Canopy_' + state + part
land_cover_layer_name = 'NLCD_' + year + '_Land_Cover_' + state + part

extent = None

################################################## Step 1: Reclass tree canopy to one landcover type 41, 42, or 43 #########################################

# Define input and output paths
output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Trees-Combined.tiff'

# Retrieve input layers from the project using iface for land cover layer
layer1 = QgsProject.instance().mapLayersByName(tree_canopy_layer_name)[0]
# Set the layers as the active layer
iface.setActiveLayer(layer1)
# Retrieve input layers from the project using iface
input_layer_tree_canopy = iface.activeLayer()

expression = (
'(A > 0) * 41 + '
'(A <= 0) * A'
)
command = [
'gdal_calc.py',
'-A', input_layer_tree_canopy.source(),
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)
result_tree_canopy_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Trees-Combined')
QgsProject.instance().addMapLayer(result_tree_canopy_layer)

print("Step 1: Calculate result_tree_canopy completed. (NLCD_" + year + "_" + state + part + "_Trees-Combined)")

############################################################# Step 2: Warp tree canopy to 4326 #############################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Tree_Canopy_4326.tiff'

layer1 = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Trees-Combined')[0]
iface.setActiveLayer(layer1)
input_layer_tree_canopy = iface.activeLayer()

processing.run(
"gdal:warpreproject",
{
'INPUT':input_layer_tree_canopy.source(),
'SOURCE_CRS':None,
'TARGET_CRS':QgsCoordinateReferenceSystem('EPSG:4326'),
'RESAMPLING':0,
'NODATA':None,
'TARGET_RESOLUTION':None,
'OPTIONS':'',
'DATA_TYPE':1,
'TARGET_EXTENT':extent,
'TARGET_EXTENT_CRS':None,
'MULTITHREADING':False,
'EXTRA':'',
'OUTPUT':output_path
}
)

result_tree_canopy_warped_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Tree_Canopy_4326')
QgsProject.instance().addMapLayer(result_tree_canopy_warped_layer)

print("Step 2: Warp tree_canopy completed. (NLCD_" + year + "_" + state + part + "_Tree_Canopy_4326)")

################################################################### Step 3: Warp land cover 4326 ############################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Land_Cover_4326.tiff'

layer1 = QgsProject.instance().mapLayersByName(land_cover_layer_name)[0]
iface.setActiveLayer(layer1)
input_layer_land_cover = iface.activeLayer()

processing.run(
"gdal:warpreproject",
{
'INPUT':input_layer_land_cover.source(),
'SOURCE_CRS':None,
'TARGET_CRS':QgsCoordinateReferenceSystem('EPSG:4326'),
'RESAMPLING':0,
'NODATA':None,
'TARGET_RESOLUTION':None,
'OPTIONS':'',
'DATA_TYPE':1,
'TARGET_EXTENT':extent,
'TARGET_EXTENT_CRS':None,
'MULTITHREADING':False,
'EXTRA':'',
'OUTPUT':output_path
}
)

result_land_cover_warped_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Land_Cover_4326')
QgsProject.instance().addMapLayer(result_land_cover_warped_layer)

print("Step 3: Warp land_cover completed. (NLCD_" + year + "_" + state + part + "_Land_Cover_4326)")

################################################################ Step 4: Combine land cover 4326 and tree canopy 4326######################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Combined_4326.tiff'

warped_tree_canopy_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Tree_Canopy_4326')[0]
warped_land_cover_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Land_Cover_4326')[0]
iface.setActiveLayer(warped_tree_canopy_layer)
input_layer_tree_canopy = iface.activeLayer()
iface.setActiveLayer(warped_land_cover_layer)
input_layer_land_cover = iface.activeLayer()

expression = (
'((B > 0) & (B != 41) & (B != 42) & (B != 43) & (A > 0)) * A + '
'((B == 41) | (B == 42) | (B == 43) | (A <= 0)) * B'
)
command = [
'gdal_calc.py',
'-A', input_layer_tree_canopy.source(),
'-B', input_layer_land_cover.source(),
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)

result_tree_canopy_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Combined_4326')
QgsProject.instance().addMapLayer(result_tree_canopy_layer)

print("Step 4: Combine tree canopy and land cover completed. (NLCD_" + year + "_" + state + part + "_Combined_4326)")

############################################################## Step 5: Replace urban and clutter with grass ###################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Grass-Only_4326.tiff'

combined_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Combined_4326')[0]
iface.setActiveLayer(combined_layer)
input_combined_layer = iface.activeLayer()

expression = (
'(A == 11) * 41 + '
'(A == 12) * 34 + '
'(A == 21) * 26 + '
'(A == 22) * 26 + '
'(A == 23) * 26 + '
'(A == 24) * 26 + '
'(A == 31) * 27 + '
'(A == 41) * 23 + '
'(A == 42) * 24 + '
'(A == 43) * 25 + '
'(A == 51) * 30 + '
'(A == 52) * 29 + '
'(A == 71) * 26 + '
'(A == 72) * 32 + '
'(A == 73) * 31 + '
'(A == 74) * 31 + '
'(A == 75) * 32 + '
'(A == 81) * 18 + '
'(A == 82) * 19 + '
'(A == 90) * 25 + '
'(A == 95) * 35'
)
command = [
'gdal_calc.py',
'-A', input_combined_layer.source(),
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)

result_grass_only_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Grass-Only_4326')
QgsProject.instance().addMapLayer(result_grass_only_layer)

print("Step 5: Calculate result_grass_only completed. (NLCD_" + year + "_" + state + part + "_Grass-Only_4326)")

################################################################## Step 6: Reclass urban 21, 22, 23 or 24 ###########################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Urban_4326.tiff'

combined_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Combined_4326')[0]
iface.setActiveLayer(combined_layer)
input_combined_layer = iface.activeLayer()

expression = (
'(A == 21)*10 + '
'(A == 22)*1 + '
'(A == 23)*1 + '
'(A == 24)*2'
)
command = [
'gdal_calc.py',
'-A', input_combined_layer.source(),
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)

result_urban_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Urban_4326')
QgsProject.instance().addMapLayer(result_urban_layer)

print("Step 6: Calculate result_urban completed. (NLCD_" + year + "_" + state + part + "_Urban_4326)")

################################################################ Step 7: Remove clutter and roads from urban ###########################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Urban-Only_4326.tiff'

layer1 = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Urban_4326')[0]
iface.setActiveLayer(layer1)
input_layer_land_cover = iface.activeLayer()

processing.run(
"grass7:r.neighbors",
{
'input':input_layer_land_cover.source(),
'selection':None,
'method':1,
'size':7,
'gauss':None,
'quantile':'',
'-c':False,
'-a':False,
'weight':'',
'output':output_path,
'GRASS_REGION_PARAMETER':extent,
'GRASS_REGION_CELLSIZE_PARAMETER':0,
'GRASS_RASTER_FORMAT_OPT':'',
'GRASS_RASTER_FORMAT_META':''
}
)

result_land_cover_warped_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Urban-Only_4326')
QgsProject.instance().addMapLayer(result_land_cover_warped_layer)

print("Step 7: Remove clutter and roads from urban. (NLCD_" + year + "_" + state + part + "_Urban-Only_4326)")

############################################################## Step 8: Combine grass only and clean urban #########################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Combined-Clean_4326.tiff'

grass_only_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Grass-Only_4326')[0]
urban_only_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Urban-Only_4326')[0]
iface.setActiveLayer(urban_only_layer)
input_urban_only = iface.activeLayer()
iface.setActiveLayer(grass_only_layer)
input_grass_only = iface.activeLayer()

expression = (
'((A > 0) & (B > 0)) * A + '
'((A <= 0) | (B <= 0)) * B + (B == 0) * 44'
)
command = [
'gdal_calc.py',
'-A', input_urban_only.source(),
'-B', input_grass_only.source(),
'--type=Byte',
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)

# Resample the combined raster to a different resolution
output_path_resampled = output_path.replace("_4326.tiff", "_HD_4326.tiff")
gdal.Warp(
output_path_resampled,
output_path,
xRes=0.00005539416,
yRes=0.00005539416
)

urban_grass_combined = QgsRasterLayer(output_path_resampled, 'NLCD_' + year + '_' + state + part + '_Combined-Clean_HD_4326')
QgsProject.instance().addMapLayer(urban_grass_combined)

print("Step 8: Combined and clean completed. (NLCD_" + year + "_" + state + part + "_Combined-Clean_HD_4326)")

################################################################## Step 9: Smooth all hi res features #################################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_HD_4326.tiff'

layer1 = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Combined-Clean_HD_4326')[0]
iface.setActiveLayer(layer1)
input_layer_land_cover = iface.activeLayer()

processing.run(
"grass7:r.neighbors",
{
'input':input_layer_land_cover.source(),
'selection':None,
'method':1,
'size':7,
'gauss':None,
'quantile':'',
'-c':False,
'-a':False,
'weight':'',
'output':output_path,
'GRASS_REGION_PARAMETER':extent,
'GRASS_REGION_CELLSIZE_PARAMETER':0,
'GRASS_RASTER_FORMAT_OPT':'',
'GRASS_RASTER_FORMAT_META':''
}
)

# Convert output to Byte using gdal_translate
converted_output_path = output_path.replace("_HD_4326.tiff", "_Final-HD_4326.tiff")
command = [
'gdal_translate',
'-ot', 'Byte', # Set output data type to Byte
output_path,
converted_output_path
]
subprocess.run(command, shell=True)

result_land_cover_warped_layer = QgsRasterLayer(converted_output_path, 'NLCD_' + year + '_' + state + part + '_Final-HD_4326')
QgsProject.instance().addMapLayer(result_land_cover_warped_layer)

print("Step 9: Smooth all hi res features completed. (NLCD_" + year + "_" + state + part + "_Final-HD_4326)")
</syntaxhighlight>

===Generating the Terrain using osgdem===
Instead of using genVPB.py, you may wish to run osgdem directly.

In the Windows/Docker platform you can send the generate tile command directly to osgdem.exe, one tile at a time.

Using the NLCD raster processing convention from above, following is the the final step after creating the raster and entering bash shell with the windows version of "docker run..."

osgdem --TERRAIN --image-ext png --RGBA --no-interpolate-imagery --disable-error-diffusion --geocentric --no-mip-mapping -t ./data/California-Southern_4326-84-hd-corrected.tiff -d ./SRTM-3/N32W115.hgt -b -115 32 -114 33 --PagedLOD -l 7 --radius-to-max-visible-distance-ratio 3 -o ./output/vpb/w120n30/w115n32/ws_w115n32.osgb

Note: the --image-ext png --RGBA flags are critical to successfully building correctly placed landclasses in the final VPB generated scenery.

If you prefer to run the scenery generation manually, running the VPB osgdem process is described in more detail here: [[Virtual Planet Builder#Running VPB]].

After doing this you should have an output directory containing files of the form <code>output/vpb/w010n50/w004n50/ws_w004n50.osgb</code>, plus a host of sub-directories. Each one of these is a 1x1 tile of terrain.

to leave the container simply type <code>exit</code>.

===Packaging the Scenery===
Once you have the terrain and line features they should be packaged in a scenery directory in vpb and Terrain sub-directories respectively. E.g.
MyCoolScenery/Terrain
MyCoolScenery/vpb
It is good practise to document the data sources used in scenery generation. Some source licenses require attribution of the original data source for anything derived, published or distributed.

To assist in fulfilling these license obligations, you can create a source.xml file in the scenery directory which includes attribution information. This will then be available from within the simulator under Help->Scenery Sources, and <u>may</u> fulfil the attribution requirements of your license. '''Note that you are responsible for fulfilling any license requirements from the data, not FlightGear'''.

The format of the file is straightforward:
<?xml version="1.0"?>
<PropertyList>
<nowiki><source></nowiki>
<name>Corine Land Cover (CLC) 2018, Version 2020_20u1</name>
<nowiki><link></nowiki><nowiki>http://web.archive.org/web/20221112175615/https://land.copernicus.eu/pan-european/corine-land-cover/clc2018?tab=metadata%2A</nowiki><nowiki></link></nowiki>
<license>GMES Open License</license>
<nowiki></source></nowiki>
<nowiki><source></nowiki>
<name>NASADEM Merged DEM Global 1 arc second V001</name>
<nowiki><link></nowiki><nowiki>https://www.earthdata.nasa.gov/</nowiki><nowiki></link></nowiki>
<license>Public Domain</license>
<nowiki></source></nowiki>
<nowiki><source></nowiki>
<name>OpenStreetMap</name>
<nowiki><link></nowiki><nowiki>https://www.openstreetmap.org/copyright</nowiki><nowiki></link></nowiki>
<license>Open Data Commons Open Database License</license>
<nowiki></source></nowiki>
</PropertyList>

Howto:Create WS3.0 terrain

2024-08-01T03:14:34Z

Wlbragg: /* Using the Python Console */

{{WS30 Navbar}}
This article provides instructions on how to generate basic WS3.0 terrain.

WS3.0 terrain consists of three parts:

# A terrain mesh consisting of a landclass texture draped over an elevation model.
# A high resolution water raster used to show water features such as rivers, lakes and coastline with more definition
# Line features such as roads and railways.

The terrain is generated by a set of tools that are packaged in a docker image for convenience.[[File:Diagram-export-21-12-2023-16 29 37.png|thumb|Basic WS3.0 Scenery Generation Process]]

== Prerequisites ==

=== Set up a Workspace ===
Create a directory with the following sub-directories:

* <code>cache</code>
* <code>data</code>
* <code>output/vpb</code>
* <code>output/Terrain</code>

=== Docker ===

# Install [https://docs.docker.com/get-started/ Docker] on your platform.
#Pull the docker image by running the following command

docker pull flightgear/ws30-vpb-generator:latest
Optionally, if you are hitting rate limits:
#Create an account on https://hub.docker.com/. (Note that you will need to click on an email verification link before you can log in for the first time)
#Run <code>docker login</code> before the '''docker pull''' command above

== Getting the base data ==
You need two pieces of data for the area of scenery you are generating:

# An elevation model (aka DEM). This indicates what altitude each point of the surface is.
# Landclass data showing what type of terrain is at each point of the surface. This is often either a Raster (effectively a texture), or vector data.

=== Elevation Model ===
Download the NASADEM elevation model for the area of scenery you wish to generate. This is available in 1x1 degree blocks from [https://lpdaac.usgs.gov/products/nasadem_hgtv001/ here], and with an interactive browser [https://search.earthdata.nasa.gov/search here].

Unzip the files into the <code>data</code> directory.

=== Landclass Raster ===
Download an landclass raster for the area of scenery you wish to generate.

* For Europe, use of [https://land.copernicus.eu/pan-european/corine-land-cover/clc2018 CORINE] is recommended.
* For the USA [https://www.mrlc.gov/viewer/ NLCD] is recommended
* Sentinel-2 data is available for the entire world via [https://livingatlas.arcgis.com/landcoverexplorer/ ESRI], but has limited set of landclasses.

Put these into the <code>data</code> directory.

More detailed terrain can be created by modifying the landclass raster, and/or generating a new raster from vector data. These processes are discussed below.

== Generating Terrain ==
To generate terrain you need to run the tools within the docker container we installed above. The docker image is like a small, independent virtual computing environment running within your system. This particular docker image has all the scenery generation tools already installed.

=== Running the docker container ===
Firstly, get the container running from the directory containing your <code>cache</code>, <code>data</code> and <code>output</code> directories:
docker run --rm --mount "type=bind,source=`pwd`/data,target=/home/flightgear/data,readonly" --mount "type=bind,source=`pwd`/output,target=/home/flightgear/output" --mount "type=bind,source=`pwd`/cache,target=/home/flightgear/cache" -it flightgear/ws30-vpb-generator:latest /bin/bash

You should now find yourself in a bash shell within your container. You should see data and output directories which are linked to the directories you created earlier:
flightgear@ddcac77f7d5e:~$ ls
cache data output bin scripts

=== Building the terrain ===
To build the terrain mesh, use the <code>genVPB.py</code> tool from inside the docker container.
Usage: genVPB.py <min-lon> <min-lat> <max-lon> <max-lat> <input-raster> [--reclass <reclass>] [--coastline <coastline>] [--hgt-dir <hgt-dir>] [--output-dir <output-dir>]
<min-lon> <min-lat> <max-lon> <max-lat> - bounding box of scenery to generated
<input-raster> - Input landclass raster
[coastline] - Optional coastline polygon data (.osm) to clip against. E.g. from <nowiki>https://osmdata.openstreetmap.de/data/land-polygons.html</nowiki>
[reclass] - Optional file containing a set of rules for reclassifying the raster similar to r.reclass. See fgmeta/ws30/mappings/
[hgt-dir] - Optional directory containing unzipped NASADEM HGT files. Defaults to '/home/flightgear/data'.
[output-dir] - Optional directory for output scenery. Defaults to '/home/flightgear/output/vpb'

For example, to generate a piece of terrain around Edinburgh (latitude 55.5, longitude 3 degrees West)
./scripts/genVPB.py -3 55 -2 56 ./data/uk_wgs84_10m_N54.tif
If you are using anything other than a CORINE raster you will need to reclassify the data to match the landclasses used by FlightGear. Those classes are defined in [https://sourceforge.net/p/flightgear/fgdata/ci/next/tree/Materials/base/landclass-mapping.xml Materials/base/landclass-mapping.xml]. You can reclassify them using the files in the [https://sourceforge.net/p/flightgear/fgmeta/ci/next/tree/ws30/mappings/ scripts/mappings] directory. E.g. to reclassify NLCD2019 data you can use <code>--reclassify ./scripts/mappings/nlcd2019.txt</code>.

genVPB.py will output the data to the output/vpb directory, in which you should find a series of files and directories.

=== Adding water ===
The terrain mesh does not have highly detailed water features - as typically the source data has a resolution of 10-25m. Water features are generated from OpenStreetMap data. To generate water features simply run the genwaterraster.py command.
Usage: ./scripts/genwaterraster.py <scenery_dir> <lon1> <lat1> <lon2> <lat2> [Cache]
<scenery_dir> Scenery directory to write to
<lon1> <lat1> Bottom left lon/lat of bounding box
<lon2> <lat2> Top right lon/lat of bounding box
[Cache] Optional directory to cache OpenStreetMap data to reduce API requests
The scenery directory should be your output/vpb directory created above. For example, to generate water for the area around Edinburgh:
./scripts/genwaterraster.py ./output/vpb -4 55 -3 56 ./cache
Inside the ./output/vpb directory there should be a set of directories and .png files.

=== Adding roads and railways ===
The terrain mesh does not have any line features - things like roads. These are generated separately from OpenStreetMap data. To generate line features simply run the genroads.py command:
Usage: ./scripts/genroads.py <scenery_dir> <lon1> <lat1> <lon2> <lat2> [Cache]
<scenery_dir> Scenery directory to write to
<lon1> <lat1> Bottom left lon/lat of bounding box
<lon2> <lat2> Top right lon/lat of bounding box
[Cache] Optional directory to cache OpenStreetMap data to reduce API requests

The scenery directory should be your output/Terrain directory created above. For example, to generate roads for the area around Edinburgh:
./scripts/genroads.py ./output/Terrain -4 55 -3 56 ./cache
Inside the ./output/Terrain directory there should be a set of directories and, .STG files text files.

==Running FlightGear with the new WS3.0 Terrain==
To test the new terrain, simply include the output directory in your scenery path and run FlightGear with the <code>--prop:/scenery/use-vpb=true</code> to enable WS3.0.

== Advanced Techniques ==
The following sections describe more complex techniques to generate higher quality WS3.0 terrain. Almost all of them involve using different data sources to generate a more detailed landclass raster before running the final scenery generation processes described above. Generating a highly detailed landclass raster is where the magic happens.

Most techniques use gdal or grass to modify the raster/vector data, typically using the QGIS program.

=== Using a different elevation model ===
If you are using another elevation model other than NASAEM, then you may need to re-project it using QGIS/gdalwarp to the WGS84 CRS (aka EPSG:4326).

=== Landclass Data Requirements ===
For any landclass data we need to ensure the data is in the correct format. That means:

# Is a Raster (geotiff) rather than Vector data. This raster will become the texture on the terrain that the terrain shaders do their magic on.
# Uses the WGS84 Coordinate Reference System. The ensures that the terrain generation step is efficient.
# Has the correct landclass values for each terrain type. We use a set of values based on the CORINE raster set, defined in [https://sourceforge.net/p/flightgear/fgdata/ci/next/tree/Materials/base/landclass-mapping.xml Materials/base/landclass-mapping.xml].

Below is a quick table showing what steps you need to take for common landclass data sources.
{| class="wikitable"
|+
!Landclass Data
!Warp to WGS84 required?
!Landclass re-classification Required?
!Raster Simplification Required?
!Conversion to Raster Required?
|-
|CORINE Raster
|Yes
|No
|No
|No
|-
|CORINE Vector
|Yes
|Yes
|No
|Yes
|-
|NLCD
|No
|Yes
|Yes
|No
|-
|Sentinel-2
|Yes
|Yes
|No
|No
|}
Conversion to Raster must be done manually. Converting to WGS84 and the correct landclasses ''can'' be done by the genVPB.py script, but slows down scenery generation. Therefore if you are planning to generate scenery multiple times it is best to pre-process the files yourself.

The easiest way to do these operations is using QGIS, which is available for most platforms. If you are scripting a toolchain, the QGIS tools include command-line equivalents for all commands.

When using QGIS, set the Project CRS to WGS84 (aka EPSG:4326). You can then add layers of Raster or Vector data from files from the <code>Layer->Add Layer</code> menu. When performing any operations, <u>always</u> write out the data to a real file so you can go back to it later. Disk space is cheap :).

=== Warping Raster Layers to WGS84 ===
genVPB.py will automatically convert to WGS84 if it detects a landclass raster with a different Coordinate System. However, you may wish to do so manually instead so all your data is on the same coordinate system.

In QGIS you can warp a raster layer to a different CRS using the Raster->Projections-Warp (Reproject) tool.

Select the following options in the dialog:

* Input Layer - Check you have selected the correct layer. The CRS is shown at the right.
* Target CRS - set to <code>EPSG:4326 - WGS84</code>, which should be your project CRS.
* Resampling Method to Use - Nearest Neighbour. (Landclass data is not like normal images. You don't want to interpolate between values.)
* Nodata value for output bands - 0.0 (This means that any data at the edges will be Ocean, usually a reasonable default)
* Advanced Parameters - No Compression
* Output data type - Byte
* Reprojected - Choose a target filename

Alternatively you can do this step from the commandline.
gdalwarp -t_srs EPSG:4326 -dstnodata 0.0 -r near -ot Byte -of GTiff -co COMPRESS=NONE -co BIGTIFF=IF_NEEDED /home/stuart/FlightGear/VPB/data/CORINE/u2018_clc2018_v2020_20u1_raster100m/DATA/U2018_CLC2018_V2020_20u1.tif /home/stuart/FlightGear/VPB/data/scratch/corine_WGS84.tif

=== Warping Vector Layers to WGS84 ===
genVPB.py will automatically convert to WGS84 if it detects a landclass raster with a different Coordinate System. However, you may wish to do so manually instead so all your data is on the same coordinate system.

In QGIS you can warp a vector layer using Vector->Data Management Tools->Reproject Layer.

Set the following options in the dialog:

* Input Layer - Check you have selected the correct layer. The CRS is shown at the right.
* Target CRS - set to <code>EPSG:4326 - WGS84</code>, which should be your project CRS.
* Reprojected - Choose a target filename

=== Reclassifying Vector Layers ===
For CORINE vector data in particular, the attributes used in the vector data are not the same as those used by the CORINE Raster data. So we need to create a new attribute on the data.
[[File:Field Calculator.png|thumb|QGIS Field Calculator]]
To do this

* Select <code>Layer->Open Attribute Table</code>. You should see a table with multiple columns. Each row represents a feature in the data.
* Click on <code>Toggle Editing Mode</code> (Ctrl + E), on the top left of the dialog
* Click on <code>New Field</code> (Ctrl+W) and create a new field called Landclass of type "Whole Number (Integer)". This will create a new column which we will populate with the correct landclass data
* Click on the <code>Open Field Calculator</code> button (Ctrl + I). (If you get an error about only being able to create Virtual fields, go back to the Layer menu, export it and open the exported file).
* Select the following options:
** Update Existing Field
** Select the Landclass field you just created.
** Copy the contents of https://sourceforge.net/p/flightgear/fgmeta/ci/next/tree/ws30/mappings/corine_vector.txt into the Expression box (without the comment lines starting with <code>#</code>). This is just some simple code to set the attribute correctly. The code should be correct for CORINE vector data. If your data is from other sources you will need to work out how you want to map your source data landclasses to the CORINE ones. [https://sourceforge.net/p/flightgear/fgdata/ci/next/tree/Materials/base/landclass-mapping.xml Materials/base/landclass-mapping.xml] can be used as a guide.
* Select <code>OK</code>. You should see that your landclass column is now populated with the landclass data.
* Select <code>Layer->Save Layer Edits</code> to save you changes

=== Creating a Raster from a Vector Layer ===
To create a Raster from a Vector Layer select <code>Raster->Conversion->Rasterize (Vector to Raster)</code>.
[[File:QGIS Rasterize (Vector to Raster).png|thumb|Creating a Raster from a Vector Layer - QGIS Rasterize]]
Select the following options:

* Input Layer - correct layer, check CRS
* Field to use for burn-in value - select the Landclass column you created above.
* Output raster size units. This is going to set the resolution of your raster. You can work out the resolution in two different ways:
** Select "Georeferenced units" and determine how many degrees each pixel is in latitude and longitude.
** Select "Pixels" and determine the size of raster you want in pixels. [https://www.nhc.noaa.gov/gccalc.shtml This] is a good calculator to help.
* Width/Horizontal Resolution. Enter the values you've calculated for the horizontal resolution (longitudinal), or the width of the raster
* Height/Vertical Resolution. Enter the values you've calculated for the vertical resolution (latitude or the height of the raster
* Output extent - Select an option from the box on the right. You can edit the text afterwards. Best practise is to create long thin strips of 1 degree latitude in height, as this makes subsequent processing much easier.
* Assign a specific nodata value to output bands - Select 0.0 for Ocean. CORINE vector data in particular has a lot of nodata for Oceans
* Advanced Parameters - No Compression
* Output data type - Byte
* Rasterized - Select a new filename

=== Simplifying a Raster Layer ===
Some Raster Landclass data (NLCD included) has too much noise - in particular large US highway systems are identified as Urban areas.

To smooth it out we can use the GRASS <code>n.neighbors</code> function from the Processing Toolbox in QGIS.

Select the following options:

* Input Layer - correct layer, check CRS
* Neighborhood operation - median. (This is not a normal image, so using an average will result in weird values)
* Neighborhood size - 5.
* Neighbors - Select a new filename.

=== Clipping a Raster Layer with OSM Data for Land (Corine) ===
The Corine dataset does not match OSM coastlines exactly. The following multi-stage process makes sure, that no Corine land-use is in the water as defined by OSM.

==== Download OSM Land Data ====

Download land polygons based on OSM data as a Shapefile from [https://osmdata.openstreetmap.de/data/land-polygons.html Land Polygons] and make sure to pick the WGS84 projected download with split polygons ("Large polygons are split, use for larger scales"). Once downloaded unzip the content into a directory.

==== Reclassifying the OSM Land Data Vector Layer ====
To do this:

* Select <code>Layer->Open Attribute Table</code>. You should see a table with multiple columns. Each row represents a feature in the data.
* Click on <code>Toggle Editing Mode</code> (Ctrl + E), on the top left of the dialogue
* Click on <code>New Field</code> (Ctrl+W) and create a new field called Landclass of type <code>Integer (32 bit)</code>. This will create a new column which we will populate with the correct land class data
* On top of the table on the left side choose "Landclass" in the drop-down menu, then input <code>2</code> into the field to the right and then press button "Update" all to the left of this field.
* Wait a bit and the close the dialogue.
* Select <code>Layer->Save Layer Edits</code> to save you changes

==== Convert the Land Data from Vector to Raster ====
Do the same as in chapter "Creating a Raster from a Vector Layer" above. The only difference is that the Input layer will be the land data polygons and you need to choose a different file name for the "Rasterized" (e.g. osm_land_scotland_5m.tif)

==== Remove Novalue Entries in the Land Data Raster ====
To do this:

* Select <code>Processing->Toolbox</code>. You should see a new box on the right side.
* Write "gdal_calc" in the search box and you should see an entry "Raster calculator". Double click on it and you will get a new dialogue window.
* In this dialogue:
** For "Input layer A" choose the raster from the previous chapter ((e.g. osm_land_scotland_5m.tif)
** In field "Calculation in gdalnumeric ..." write: <code>greater(A,0) * A</code>
** In field "Output raster type" choose <code>Byte</code>
** In field "Advanced Parameters" choose Profile <code>No compression</code>
** In field "Additional command-line parameters" write: <code>--hideNoData</code>
** In field "Calculated" choose a file (e.g. osm_land_scotland_allvalues_5m.tif)
** (In the "GDAL/OGR console call" it will have something similar to the follwing - just with different paths: <code>gdal_calc.py --overwrite --calc "greater(A ,0) * A" --format GTiff --type Byte -A /home/vanosten/custom-fg-scenery/data/osm_land_scotland_5m.tif --A_band 1 --co COMPRESS=NONE --co BIGTIFF=IF_NEEDED --hideNoData --outfile /home/vanosten/custom-fg-scenery/data/osm_land_scotland_allvalues_5m.tif</code>
** Press the "Run" button - and when complete close the dialogue.

You should now see a map only black and white. You can check for correctness by pressing <code>CTRL+SHIFT+I</code> to get a cursor with an arrow and an "i". First make sure the new raster is selected on the left side. Next click on the sea/ocean and then check in the "Identify Results" window on the right that the value is < 2. The click on the land and check that the value is 2.

==== Create the Final Clipped Corine Raster Against OSM Land Data =====
Do the following:

* In QGIS make sure that you have only the following two layers: the basis Corine raster (see chapter "Creating a Raster from a Vector Layer" - here e.g. corine_raster_scotland_5m.tif) and plus the raster from the previous step (e.g. osm_land_scotland_allvalues_5m.tif)
* Select <code>Raster->Raster Calculator ...</code> and a corresponding dialogue will open showing on the left hand side the two rasters.
* Choose a new "Output layer" (e.g. corine_raster_scotland_clipped_5m.tif).
* In the "Raster Calculator Expression" field input: <code>if ("osm_land_scotland_all_data_5m@1" < 2, 44, "corine_raster_scotland_5m@1")</code>
* Press button "OK" and wait a while (you will see a new dialogue with showing the progress.

Done. You have now a raster (e.g. corine_raster_scotland_clipped_5m.tif) which does not have land in areas, where OSM data has sea/ocean.

=== Reclassifying a Raster Layer ===
WS3.0 uses CORINE landclass values. If using data from other sources it needs to be reclassified to the correct values. genVPB.py has an option to do this, but you may wish to do so manually.

To do this select <code>GRASS->Raster->r.reclass</code> from the Processing Toolbox.

Select the following options:

* Input Raster Layer - correct layer, check CRS
* Reclass rules text - copy in the contents of https://sourceforge.net/p/flightgear/fgmeta/ci/next/tree/ws30/mappings/nlcd2019.txt. Or an appropriate mapping from your landclass data to CORINE. Note that you can also reference a file using the "File containing reclass rules" option. Note a mapping of 22 24 = 1 is the same as 22 and 24 = 1. For a range of 22 to 24 use 22 23 24 = 1.
* Reclassified - Select a new filename.

(If this doesn't work a similar function is available in the Processing Toolbox under <code>Raster analysis->Reclassify by table</code>. However this doesn't save your table once you close the dialog, and entries have to be manually entered individually which takes a lot of effort)

=== Step By Step Procedure for Processing NLCD for the USA using the Raster Calculator, Up sampling and GRASS r.neighbors===

We will use a predetermined file naming convention throughout this procedure for simplicity. You can use a naming convention of your choice.

Obtain your NLCD_2019_Land_Cover from source of choice, https://www.mrlc.gov/viewer/

NLCD_2019_Land_Cover data = NLCD_2019_Land_Cover_California-Southern.tiff

Obtain your NLCD_2016_Tree_Canopy from source of choice, https://www.mrlc.gov/viewer/

NLCD_2016_Tree_Canopy data = NLCD_2016_Tree_Canopy_California-Southern.tiff

Make Deciduous Layer 255 to 0 EPSG:5070 - NAD83 - Raster Calculator

("NLCD_2016_Tree_Canopy_California-Southern@1" = 255) * 0 + ("NLCD_2016_Tree_Canopy_California-Southern@1" != 255) * "NLCD_2016_Tree_Canopy_California-Southern@1"

Output CRS = EPSG:5070 - NAD83

Output Layer = California-Southern_deciduous-coast-clipped.tiff

OK

Make Tree Canopy Layer EPSG:5070 - NAD83 - Raster Calculator

Use the prevailing tree type, for Deciduous use * 41, for Evergreen use * 42, for MixedForest use 43.

Later we add this Tree Canopy layer to the raster where there is currently no trees in the NLCD

If ocean bordering with index 255 -

("California-Southern_deciduous-coast-clipped@1" != 0) * 42

Else

("NLCD_2016_Tree_Canopy_California-Southern@1" != 0) * 42

Output CRS = EPSG:5070 - NAD83

Output Layer = California-Southern_deciduous.tiff

OK

Make Clipped Ocean Frontage EPSG:5070 - NAD83 - GIMP

Source = NLCD_2019_Land_Cover_California-Southern.tiff

Flood fill index 11 (blue) with 0 (black)

Export tiff = California-Southern_coast-clipped.tiff

Warp NLCD Deciduous and Base NLCD from EPSG:5070 - NAD83 / Conus Albers NLCD to Project CRS: EPSG:4326 - WGS 84 - WARP

If ocean bordering with index 255

Warp California-Southern_coast-clipped.tiff

Else

Warp NLCD_2016_Tree_Canopy_California-Southern.tiff

and

Warp Deciduous NLCD from EPSG:5070 - NAD83 / Conus Albers NLCD to Project CRS: EPSG:4326 - WGS 84

Warp California-Southern_deciduous.tiff

Raster -> Projection - > Warp

Target CRS = EPSG:4326 - WGS 84

Resampling method to use = Nearest Neighbour

Output data type = byte

Reprojected = California-Southern_4326-84.tiff and California-Southern_deciduous_4326-84.tiff

Run

Combine "California-Southern_deciduous_4326-84.@1" and "California-Southern_4326-84@1" - Raster Calculator

This step adds the Tree Canopy layer to the existing NLCD raster only where there is currently no existing tree land cover.

(("California-Southern_4326-84@1" > 0 AND "California-Southern_4326-84@1" != 41 AND "California-Southern_4326-84@1" != 42 AND "California-Southern_4326-84@1" != 43) AND "California-Southern_deciduous_4326-84@1" > 0) * "California-Southern_deciduous_4326-84@1" + ("California-Southern_4326-84@1" = 41 OR "California-Southern_4326-84@1" = 42 OR "California-Southern_4326-84@1" = 43 OR "California-Southern_deciduous_4326-84@1" <= 0) * "California-Southern_4326-84@1"

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_adjusted.tiff

OK

Re-class urban 21, 22, 23, 24 = grass 26 - Raster Calculator

This step is being used to erase all the man made clutter on the raster with grassland. Instead of grassland you can replace it with any type of prevailing groundcover in the area your processing.

For example this places a road easement of grass on all the road clutter that was originally in the raster. If your processing a desert area you might want to use dirt or sand instead of grassland. For Arctic areas maybe use tundra.

Part of the clutter gets replaced later in the process.

("California-Southern_adjusted@1" = 11) * 41 + ("California-Southern_adjusted@1" = 12) * 34 + ("California-Southern_adjusted@1" = 21) * 26 + ("California-Southern_adjusted@1" = 22) * 26 + ("California-Southern_adjusted@1" = 23) * 26 + ("California-Southern_adjusted@1" = 24) * 26 + ("California-Southern_adjusted@1" = 31) * 27 + ("California-Southern_adjusted@1" = 41) * 26 + ("California-Southern_adjusted@1" = 42) * 24 + ("California-Southern_adjusted@1" = 43) * 25 + ("California-Southern_adjusted@1" = 51) * 30 + ("California-Southern_adjusted@1" = 52) * 29 + ("California-Southern_adjusted@1" = 71) * 26 + ("California-Southern_adjusted@1" = 72) * 32 + ("California-Southern_adjusted@1" = 73) * 31 + ("California-Southern_adjusted@1" = 74) * 31 + ("California-Southern_adjusted@1" = 75) * 32 + ("California-Southern_adjusted@1" = 81) * 18 + ("California-Southern_adjusted@1" = 82) * 19 + ("California-Southern_adjusted@1" = 90) * 25 + ("California-Southern_adjusted@1" = 95) * 35

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_reclass-grass.tiff

OK

Make Urban layer

("California-Southern_adjusted@1" = 21) * 10 + ("California-Southern_adjusted@1" = 22) * 1 + ("California-Southern_adjusted@1" = 23) * 1 + ("California-Southern_adjusted@1" = 24) * 2

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_urban.tiff

OK

Remove the road clutter with r.neighbor

GRASS/Raster/r.neighbor

Neighborhood operation = median

Neighborhood size = 7

Neighbors = California-Southern_urban-only.tiff

Run

Combine "California-Southern_reclass-urban@1" and "California-Southern_reclass-grass@1"

This step replaces the "cleaned up" urban areas back into the raster.

("California-Southern_urban-only@1" < 1) * "California-Southern_reclass-grass@1" + ("California-Southern_urban-only@1" != 0) * "California-Southern_urban-only@1"

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_adjusted-combined.tiff

OK

Up sample California-Southern_adjusted-combined.tiff 2x, 4x or 8x resolution 0.00008309125

Active Layer = California-Southern_adjusted-combined.tiff

Layer/Save As

2x Horizontal = 0.00008309125

2x Vertical = 0.00008309125

4x Horizontal = 0.00010637625

4x Vertical = 0.00010637625

8x Horizontal = 0.00005318812

8x Vertical = 0.00005318812

File Name = California-Southern_final-prep-4x.tiff

OK

Simplify California-Southern_4326-84-hd.tiff with r.neighbor

GRASS/Raster/r.neighbor

Neighborhood operation = median

Neighborhood size = 7

Neighbors = California-Southern_urban-only.tiff

Run

Reclass index 0 to 44, leave the rest the same

("California-Southern_4326-84-hd@1" = 0) * 44 + ("California-Southern_4326-84-hd@1" != 0) * "California-Southern_4326-84-hd@1"

Output CRS = EPSG:4326 - WGS 84

Output Layer = California-Southern_4326-84-hd-corrected.tiff

OK

<big>Optional HD Fresh Water Option</big>
*Step 1 - Obtain and load hi resolution vector layer
Make sure that the vector layer and the raster layer you will eventually merge to have the same projection. ** Extent can be different if you use the option below.
#Use Top Menu: "Raster" -> "Conversion" -> "Rasterize (vector to raster)" or Processing Toolbox: GDAL -> "Vector conversion" -> "Rasterize (vector to raster)"
# Input layer = Southern-California_water_4326-84
#Fixed value to burn = 41 (water)
# Output raster size units = "Georeferenced units"
#Width/Horizontal resolution = 0.00008309125 (4x the base NLCD resolution )
#Height/Vertical resolution = 0.00008309125 (4x the base NLCD resolution)
#<nowiki>** Output extent = Select the raster layer you will eventually merge with as the "Output extent".</nowiki>
#Under "Rasterized" choose a new filename to save the new raster to. We'll use Southern-California_4326-84-hd-water.tiff.

Click "Run" to generate the hi resolution raster layer.

*Step 2 - Reclass hi res, smoothed water to grass
# Open Raster Calculator
#Under "Output layer" choose a new filename to save the new raster to. We'll use Southern-California_4326-84-hd-nowater.tiff.
#Copy the following into the "Raster Calculator Expression" box. Where Southern-California_4326-84-hd is the name of your hi resolution, smoothed NLCD that includes the water data.
"Southern-California_4326-84-hd@1" * ("Southern-California_4326-84-hd@1" != 41) + 26 * ("Southern-California_4326-84-hd@1 = 41")

Click on "OK". When finished you will have a new raster with the water layer changed to grassland. Another choice would be to change it to sand.

*Step 3 - Combine the hi resolution no water raster and the hi resolution water raster.
# Use Top Menu: "Raster" -> "Miscellaneous" -> "Merge" or Processing Toolbax: GDAL -> Raster Miscellaneous -> Merge
#Input layers = Select "Southern-California_4326-84-hd-nowater@1" and Southern-California_4326-84-hd-water.tiff
#Output data type = "byte".
# Under "Merged" choose a new filename to save the new raster to. We'll use Southern-California_4326-84-hdwater.tiff.

Click "Run" to generate the new merged hi resolution raster layer.

==='''Optional Python script''' to process NLCD for the USA===
You can optionally process the NLCD by loading this script in the python console of the QGIS Desktop program.

If you save the NLCD tiff's using the same consistent naming convention used in this example it is fairly simple to generate final WS3.0 tiff's to use to generate the scenery.

You will start with an original tree layer and a land cover layer from the MRLC.gov site. See the section "Obtaining NLCD" in the above "Step By Step Procedure for Processing NLCD for the USA using the Raster Calculator" topic.

For example, to start with you will have a couple files named...

NLCD_2021_Tree_Canopy_NewYork.tiff

NLCD_2021_Land_Cover_NewYork.tiff

In the QGIS program the layer names you start with will be...

NLCD_2021_Tree_Canopy_NewYork

NLCD_2021_Land_Cover_NewYork

which will refer to the tiff files.

Running this script will generate a bunch of intermediate files including the last file, that will be the final finished file for running in the VPB script or VPB build environment Docker image.

The last file for this example for New York would be, NLCD_2021_Land_Cover_NewYork_Final-HD_4326.tiff<syntaxhighlight lang="python">
from qgis.core import QgsProject, QgsRasterLayer, QgsCoordinateReferenceSystem, QgsProcessingException
from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry
from osgeo import gdal, osr
import os
import numpy
import numpy as np
import processing
import subprocess

# Define input layer names, change according to your file names
year = '2021';
state = 'Kentucky';
part = '';

tree_canopy_layer_name = 'NLCD_' + year + '_Tree_Canopy_' + state + part
land_cover_layer_name = 'NLCD_' + year + '_Land_Cover_' + state + part

extent = None

################################################## Step 1: Reclass tree canopy to one landcover type 41, 42, or 43 #########################################

# Define input and output paths
output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Trees-Combined.tiff'

# Retrieve input layers from the project using iface for land cover layer
layer1 = QgsProject.instance().mapLayersByName(tree_canopy_layer_name)[0]
# Set the layers as the active layer
iface.setActiveLayer(layer1)
# Retrieve input layers from the project using iface
input_layer_tree_canopy = iface.activeLayer()

expression = (
'(A > 0) * 41 + '
'(A <= 0) * A'
)
command = [
'gdal_calc.py',
'-A', input_layer_tree_canopy.source(),
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)
result_tree_canopy_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Trees-Combined')
QgsProject.instance().addMapLayer(result_tree_canopy_layer)

print("Step 1: Calculate result_tree_canopy completed. (NLCD_" + year + "_" + state + part + "_Trees-Combined)")

############################################################# Step 2: Warp tree canopy to 4326 #############################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Tree_Canopy_4326.tiff'

layer1 = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Trees-Combined')[0]
iface.setActiveLayer(layer1)
input_layer_tree_canopy = iface.activeLayer()

processing.run(
"gdal:warpreproject",
{
'INPUT':input_layer_tree_canopy.source(),
'SOURCE_CRS':None,
'TARGET_CRS':QgsCoordinateReferenceSystem('EPSG:4326'),
'RESAMPLING':0,
'NODATA':None,
'TARGET_RESOLUTION':None,
'OPTIONS':'',
'DATA_TYPE':1,
'TARGET_EXTENT':extent,
'TARGET_EXTENT_CRS':None,
'MULTITHREADING':False,
'EXTRA':'',
'OUTPUT':output_path
}
)

result_tree_canopy_warped_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Tree_Canopy_4326')
QgsProject.instance().addMapLayer(result_tree_canopy_warped_layer)

print("Step 2: Warp tree_canopy completed. (NLCD_" + year + "_" + state + part + "_Tree_Canopy_4326)")

################################################################### Step 3: Warp land cover 4326 ############################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Land_Cover_4326.tiff'

layer1 = QgsProject.instance().mapLayersByName(land_cover_layer_name)[0]
iface.setActiveLayer(layer1)
input_layer_land_cover = iface.activeLayer()

processing.run(
"gdal:warpreproject",
{
'INPUT':input_layer_land_cover.source(),
'SOURCE_CRS':None,
'TARGET_CRS':QgsCoordinateReferenceSystem('EPSG:4326'),
'RESAMPLING':0,
'NODATA':None,
'TARGET_RESOLUTION':None,
'OPTIONS':'',
'DATA_TYPE':1,
'TARGET_EXTENT':extent,
'TARGET_EXTENT_CRS':None,
'MULTITHREADING':False,
'EXTRA':'',
'OUTPUT':output_path
}
)

result_land_cover_warped_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Land_Cover_4326')
QgsProject.instance().addMapLayer(result_land_cover_warped_layer)

print("Step 3: Warp land_cover completed. (NLCD_" + year + "_" + state + part + "_Land_Cover_4326)")

################################################################ Step 4: Combine land cover 4326 and tree canopy 4326######################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Combined_4326.tiff'

warped_tree_canopy_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Tree_Canopy_4326')[0]
warped_land_cover_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Land_Cover_4326')[0]
iface.setActiveLayer(warped_tree_canopy_layer)
input_layer_tree_canopy = iface.activeLayer()
iface.setActiveLayer(warped_land_cover_layer)
input_layer_land_cover = iface.activeLayer()

expression = (
'((B > 0) & (B != 41) & (B != 42) & (B != 43) & (A > 0)) * A + '
'((B == 41) | (B == 42) | (B == 43) | (A <= 0)) * B'
)
command = [
'gdal_calc.py',
'-A', input_layer_tree_canopy.source(),
'-B', input_layer_land_cover.source(),
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)

result_tree_canopy_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Combined_4326')
QgsProject.instance().addMapLayer(result_tree_canopy_layer)

print("Step 4: Combine tree canopy and land cover completed. (NLCD_" + year + "_" + state + part + "_Combined_4326)")

############################################################## Step 5: Replace urban and clutter with grass ###################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Grass-Only_4326.tiff'

combined_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Combined_4326')[0]
iface.setActiveLayer(combined_layer)
input_combined_layer = iface.activeLayer()

expression = (
'(A == 11) * 41 + '
'(A == 12) * 34 + '
'(A == 21) * 26 + '
'(A == 22) * 26 + '
'(A == 23) * 26 + '
'(A == 24) * 26 + '
'(A == 31) * 27 + '
'(A == 41) * 23 + '
'(A == 42) * 24 + '
'(A == 43) * 25 + '
'(A == 51) * 30 + '
'(A == 52) * 29 + '
'(A == 71) * 26 + '
'(A == 72) * 32 + '
'(A == 73) * 31 + '
'(A == 74) * 31 + '
'(A == 75) * 32 + '
'(A == 81) * 18 + '
'(A == 82) * 19 + '
'(A == 90) * 25 + '
'(A == 95) * 35'
)
command = [
'gdal_calc.py',
'-A', input_combined_layer.source(),
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)

result_grass_only_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Grass-Only_4326')
QgsProject.instance().addMapLayer(result_grass_only_layer)

print("Step 5: Calculate result_grass_only completed. (NLCD_" + year + "_" + state + part + "_Grass-Only_4326)")

################################################################## Step 6: Reclass urban 21, 22, 23 or 24 ###########################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Urban_4326.tiff'

combined_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Combined_4326')[0]
iface.setActiveLayer(combined_layer)
input_combined_layer = iface.activeLayer()

expression = (
'(A == 21)*10 + '
'(A == 22)*1 + '
'(A == 23)*1 + '
'(A == 24)*2'
)
command = [
'gdal_calc.py',
'-A', input_combined_layer.source(),
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)

result_urban_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Urban_4326')
QgsProject.instance().addMapLayer(result_urban_layer)

print("Step 6: Calculate result_urban completed. (NLCD_" + year + "_" + state + part + "_Urban_4326)")

################################################################ Step 7: Remove clutter and roads from urban ###########################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Urban-Only_4326.tiff'

layer1 = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Urban_4326')[0]
iface.setActiveLayer(layer1)
input_layer_land_cover = iface.activeLayer()

processing.run(
"grass7:r.neighbors",
{
'input':input_layer_land_cover.source(),
'selection':None,
'method':1,
'size':7,
'gauss':None,
'quantile':'',
'-c':False,
'-a':False,
'weight':'',
'output':output_path,
'GRASS_REGION_PARAMETER':extent,
'GRASS_REGION_CELLSIZE_PARAMETER':0,
'GRASS_RASTER_FORMAT_OPT':'',
'GRASS_RASTER_FORMAT_META':''
}
)

result_land_cover_warped_layer = QgsRasterLayer(output_path, 'NLCD_' + year + '_' + state + part + '_Urban-Only_4326')
QgsProject.instance().addMapLayer(result_land_cover_warped_layer)

print("Step 7: Remove clutter and roads from urban. (NLCD_" + year + "_" + state + part + "_Urban-Only_4326)")

############################################################## Step 8: Combine grass only and clean urban #########################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_Combined-Clean_4326.tiff'

grass_only_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Grass-Only_4326')[0]
urban_only_layer = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Urban-Only_4326')[0]
iface.setActiveLayer(urban_only_layer)
input_urban_only = iface.activeLayer()
iface.setActiveLayer(grass_only_layer)
input_grass_only = iface.activeLayer()

expression = (
'((A > 0) & (B > 0)) * A + '
'((A <= 0) | (B <= 0)) * B + (B == 0) * 44'
)
command = [
'gdal_calc.py',
'-A', input_urban_only.source(),
'-B', input_grass_only.source(),
'--type=Byte',
'--outfile', output_path,
'--calc', expression
]
subprocess.run(command, shell=True)

# Resample the combined raster to a different resolution
output_path_resampled = output_path.replace("_4326.tiff", "_HD_4326.tiff")
gdal.Warp(
output_path_resampled,
output_path,
xRes=0.00005539416,
yRes=0.00005539416
)

urban_grass_combined = QgsRasterLayer(output_path_resampled, 'NLCD_' + year + '_' + state + part + '_Combined-Clean_HD_4326')
QgsProject.instance().addMapLayer(urban_grass_combined)

print("Step 8: Combined and clean completed. (NLCD_" + year + "_" + state + part + "_Combined-Clean_HD_4326)")

################################################################## Step 9: Smooth all hi res features #################################################################

output_path = 'G:/Scenery/ws3.0/' + state + '/data/NLCD_' + year + '_' + state + part + '_HD_4326.tiff'

layer1 = QgsProject.instance().mapLayersByName('NLCD_' + year + '_' + state + part + '_Combined-Clean_HD_4326')[0]
iface.setActiveLayer(layer1)
input_layer_land_cover = iface.activeLayer()

processing.run(
"grass7:r.neighbors",
{
'input':input_layer_land_cover.source(),
'selection':None,
'method':1,
'size':7,
'gauss':None,
'quantile':'',
'-c':False,
'-a':False,
'weight':'',
'output':output_path,
'GRASS_REGION_PARAMETER':extent,
'GRASS_REGION_CELLSIZE_PARAMETER':0,
'GRASS_RASTER_FORMAT_OPT':'',
'GRASS_RASTER_FORMAT_META':''
}
)

# Convert output to Byte using gdal_translate
converted_output_path = output_path.replace("_HD_4326.tiff", "_Final-HD_4326.tiff")
command = [
'gdal_translate',
'-ot', 'Byte', # Set output data type to Byte
output_path,
converted_output_path
]
subprocess.run(command, shell=True)

result_land_cover_warped_layer = QgsRasterLayer(converted_output_path, 'NLCD_' + year + '_' + state + part + '_Final-HD_4326')
QgsProject.instance().addMapLayer(result_land_cover_warped_layer)

print("Step 9: Smooth all hi res features completed. (NLCD_" + year + "_" + state + part + "_Final-HD_4326)")
</syntaxhighlight>

===Generating the Terrain using osgdem===
Instead of using genVPB.py, you may wish to run osgdem directly.

In the Windows/Docker platform you can send the generate tile command directly to osgdem.exe, one tile at a time.

Using the NLCD raster processing convention from above, following is the the final step after creating the raster and entering bash shell with the windows version of "docker run..."

osgdem --TERRAIN --image-ext png --RGBA --no-interpolate-imagery --disable-error-diffusion --geocentric --no-mip-mapping -t ./data/California-Southern_4326-84-hd-corrected.tiff -d ./SRTM-3/N32W115.hgt -b -115 32 -114 33 --PagedLOD -l 7 --radius-to-max-visible-distance-ratio 3 -o ./output/vpb/w120n30/w115n32/ws_w115n32.osgb

Note: the --image-ext png --RGBA flags are critical to successfully building correctly placed landclasses in the final VPB generated scenery.

If you prefer to run the scenery generation manually, running the VPB osgdem process is described in more detail here: [[Virtual Planet Builder#Running VPB]].

After doing this you should have an output directory containing files of the form <code>output/vpb/w010n50/w004n50/ws_w004n50.osgb</code>, plus a host of sub-directories. Each one of these is a 1x1 tile of terrain.

to leave the container simply type <code>exit</code>.

===Packaging the Scenery===
Once you have the terrain and line features they should be packaged in a scenery directory in vpb and Terrain sub-directories respectively. E.g.
MyCoolScenery/Terrain
MyCoolScenery/vpb
It is good practise to document the data sources used in scenery generation. Some source licenses require attribution of the original data source for anything derived, published or distributed.

To assist in fulfilling these license obligations, you can create a source.xml file in the scenery directory which includes attribution information. This will then be available from within the simulator under Help->Scenery Sources, and <u>may</u> fulfil the attribution requirements of your license. '''Note that you are responsible for fulfilling any license requirements from the data, not FlightGear'''.

The format of the file is straightforward:
<?xml version="1.0"?>
<PropertyList>
<nowiki><source></nowiki>
<name>Corine Land Cover (CLC) 2018, Version 2020_20u1</name>
<nowiki><link></nowiki><nowiki>http://web.archive.org/web/20221112175615/https://land.copernicus.eu/pan-european/corine-land-cover/clc2018?tab=metadata%2A</nowiki><nowiki></link></nowiki>
<license>GMES Open License</license>
<nowiki></source></nowiki>
<nowiki><source></nowiki>
<name>NASADEM Merged DEM Global 1 arc second V001</name>
<nowiki><link></nowiki><nowiki>https://www.earthdata.nasa.gov/</nowiki><nowiki></link></nowiki>
<license>Public Domain</license>
<nowiki></source></nowiki>
<nowiki><source></nowiki>
<name>OpenStreetMap</name>
<nowiki><link></nowiki><nowiki>https://www.openstreetmap.org/copyright</nowiki><nowiki></link></nowiki>
<license>Open Data Commons Open Database License</license>
<nowiki></source></nowiki>
</PropertyList>