TransWikia.com

How to Extract (Magnetic) XYZ Data Points from Color (Airborne Magnetic Survey) Image and Legend

Mathematica Asked by Inceptionalist on May 2, 2021

I’ve been trying to work my way through this by working through search results like:

Data extraction from a picture of a graph

Numeric data from a color map

Extract bathymetry data from a map

I can’t seem to make it work.

Basically, I have a GEOTIFF of the magnetic survey flown but I don’t have the raw data. The GEOTIFF is here: Airborne Magnetic Survey

And the legend is here: Airborne Magnetic Color Legend

What I’m trying to do is to extract the Z data from the image (ie. the value in nT). I can figure out how to get X,Y but I’m struggling to figure out how to obtain the value. Any advice would be huge. Ideally, I’d like to be able to get data from the map into a CSV I can work with. Thanks!

One Answer

legend = Import["https://i.imgur.com/onzvkyP.png"];

(* extract a slice of columns containing the colours. Rotate the legend on its side. *)
legadjusted = ImageRotate[ImageCrop[ImageTake[legend, All, {40, 45}]], -Pi/2];

(* we generate a list of points that fall in between the boxes in the legend *)
coords = Table[{x, 5}, {x, 16, 398, 9.5}];

(* illustrates where I've put the points *)
Show[legadjusted, Graphics[{Point[coords]}]]

(* extract the colours at the pixel coordinates, ignore the last few bad ones*)
colours = PixelValue[legadjusted, coords][[;; -3]];

(* create a range of values *)
values = N@Range[54443, 55845, (55845 - 54443 + 1)/Length[colours]];

(* load the image and build a NearestFunction that maps colours to indices in the colour map*)
img = Import["https://i.stack.imgur.com/bi16H.jpg"];
nf = Nearest[colours -> "Index"];

(* gets the index of the nearest colours and the corresponding value *)
col2val = Function[{c}, values[[First@nf[c]]]];


(* we are done! now get a grid of values by mapping over the image *)
datagrid = Map[col2val, ImageData[img], {2}];

(* we can even turn it into monochrome if we want *)
Image[datagrid] // ImageAdjust

(* or even make a 3D landscape *)
ListPlot3D[Reverse@datagrid, PlotTheme -> "Web"]

landscape

Answered by flinty on May 2, 2021

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP