# How to create a RGB plot from three different TIF files?

**URL:** <https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001>\
**Category:** Geo\
**Tags:** plotting\
**Created:** [December 22, 2022, 1:07pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001 "2022-12-22T13:07:36Z")\
**Posts on this page:** 18\
**Page:** 1

<div class="post-metadata">

**Author:** ![EmanuelCastanho](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/emanuelcastanho/32/45245_2.png) [@EmanuelCastanho](https://discourse.julialang.org/u/EmanuelCastanho)\
**Post date:** [December 22, 2022, 1:07pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/1 "2022-12-22T13:07:36Z")

</div>

I have three TIF files that represent the Red, Green and Blue channels, they are pós-processed Sentinel-2 images.

I read them with ArchGDAL and stacked them into a Float32 (554, 649, 3) matrix.

How can I plot this stacked matrix as a RGB image to then overlay a uInt8 matrix of the same size on the same plot (preferably a heat map)?

---

<div class="post-metadata">

**Author:** ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)\
**Post date:** [December 22, 2022, 1:26pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/2 "2022-12-22T13:26:59Z")

</div>

You can construct a matrix of RGB objects from Colors.jl and use Makie commands to plot like `image`:

[https://docs.makie.org/v0.19.0/examples/plotting\_functions/image/index.html](https://docs.makie.org/v0.19.0/examples/plotting_functions/image/index.html)

In MeshViz.jl we wrap some of these low-level functions into more high-level recipes for geospatial data:

> **[GitHub - JuliaGeometry/MeshViz.jl: Makie.jl recipes for visualization of...](https://github.com/JuliaGeometry/MeshViz.jl)**
>
> Makie.jl recipes for visualization of Meshes.jl. Contribute to JuliaGeometry/MeshViz.jl development by creating an account on GitHub.

You can watch my JuliaCon talk to get an idea of what is possible:

[![](https://global.discourse-cdn.com/julialang/original/3X/8/6/86a1676a1da38b4db999abd19ef80e6819092e21.jpeg "Geostatistical Learning | Júlio Hoffimann | JuliaCon 2021") ](https://www.youtube.com/watch?v=75A6zyn5pIE)

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [December 22, 2022, 6:01pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/3 "2022-12-22T18:01:09Z")

</div>

You may also find this example useful. It creates a true color image from a Landsat8 scene.

[https://www.generic-mapping-tools.org/RemoteS.jl/dev/gallery/L8cube\_img/remotes\_L8\_cube\_img/](https://www.generic-mapping-tools.org/RemoteS.jl/dev/gallery/L8cube_img/remotes_L8_cube_img/)

---

<div class="post-metadata">

**Author:** ![tim.holy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tim.holy/32/52_2.png) [@tim.holy](https://discourse.julialang.org/u/tim.holy)\
**Post date:** [December 23, 2022, 8:42pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/4 "2022-12-23T20:42:58Z")

</div>

With `ImageCore` you can also use `colorview(RGB, redchannel, greenchannel, bluechannel)`.

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [December 24, 2022, 1:46am UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/5 "2022-12-24T01:46:47Z")

</div>

Yes, joining three bands in a RGB image is simple, but if images comes from sentinel satellite it’s expected that they retain the geo-referentiation and that’s what the `RemoteS` procedure does.

---

<div class="post-metadata">

**Author:** ![EmanuelCastanho](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/emanuelcastanho/32/45245_2.png) [@EmanuelCastanho](https://discourse.julialang.org/u/EmanuelCastanho)\
**Post date:** [December 27, 2022, 2:56pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/6 "2022-12-27T14:56:13Z")

</div>

I was checking the `RemoteS.jl` repository and following the example there, however I have the following error:  
`MethodError: no method matching gmtread(::GMTgrid{Float32, 2})`

This error only appears when using` Irgb = truecolor(Ir, Ig, Ib);`

The following lines work without a problem:

```julia
Ir = gmtread("Inputs/Bands_Indices-S2/B04.tif");
Ig = gmtread("Inputs/Bands_Indices-S2/B03.tif");
Ib = gmtread("Inputs/Bands_Indices-S2/B02.tif");

```

The typeof Ir is `2267×4249 GMTgrid{Float32, 2}`.

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [December 27, 2022, 4:54pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/7 "2022-12-27T16:54:03Z")

</div>

Hmm, I see. I think that the `truecolor` function is expecting an unsigned integer but you are passing them floats. GMT has an `imagesc` function that should be what you need to convert to that `truecolor` expects.

```julia
Ir_img = imagesc(Ir);
Ig_img = imagesc(Ig);
Ib_img = imagesc(Ib);

Irgb = truecolor(Ir_img, Ig_img, Ib_img);

```

---

<div class="post-metadata">

**Author:** ![EmanuelCastanho](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/emanuelcastanho/32/45245_2.png) [@EmanuelCastanho](https://discourse.julialang.org/u/EmanuelCastanho)\
**Post date:** [December 27, 2022, 6:30pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/8 "2022-12-27T18:30:14Z")

</div>

So, after installing `ghostscript`, I was able to obtain a plot. I was not expecting this output.

 ![Screenshot 2022-12-27 at 17.25.44](https://global.discourse-cdn.com/julialang/original/3X/e/d/edecb258f6bee3b6b24e368a949a4d968e6b21c9.jpeg)

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [December 27, 2022, 6:47pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/9 "2022-12-27T18:47:06Z")

</div>

Without knowing what’s in the original files it starts to be difficult (can you make them available somewhere?)  
Does this plot a meaningful image?

```julia
imshow(Ir)

```

---

<div class="post-metadata">

**Author:** ![EmanuelCastanho](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/emanuelcastanho/32/45245_2.png) [@EmanuelCastanho](https://discourse.julialang.org/u/EmanuelCastanho)\
**Post date:** [December 29, 2022, 10:34am UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/10 "2022-12-29T10:34:34Z")

</div>

You are right, here is an example: [Example.zip - Google Drive](https://drive.google.com/file/d/1elBnvvpL1_AEQMv8vf5MblYRouhBJJgD/view?usp=sharing)

The problem seems to be when using `truecolor`.  
The outputs of `gmtread` and `imagesc` for Ir, looks right (Ir\_img is flipped):

 ![Ir](https://global.discourse-cdn.com/julialang/original/3X/1/6/16677206f0d64bfcee4554a6ed04bae16e4326c2.png)  
 ![Ir_img](https://global.discourse-cdn.com/julialang/original/3X/e/8/e8e4b3899af8fab13ac36ac655b3a9b44c6223ec.png)

Also, I don’t know if it is better to open an issue on the repository and only post the solution here.

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [December 29, 2022, 2:54pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/11 "2022-12-29T14:54:59Z")

</div>

Yes, it’s better to continue this over an issue on RemoteS so please open one. But notice that your float grids have very little dynamic

```julia
v_min: 0.039468023926 v_max: 0.159227579832

```

and I’m guessing that the bright dots are due to clouds. Did you do a processing like the one explained in this [tutorial](https://www.generic-mapping-tools.org/GMTjl_doc/tutorials/Landsat8/histogram_stretch/)? The histogram stretch step is an absolute must.

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [December 29, 2022, 4:34pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/12 "2022-12-29T16:34:47Z")

</div>

It turned out those were simple fixes (after finding the reason). The issue was that GMT accepts grids and images with several different memory layouts (row & column major, top-down or bottom-up, band-line-pixel interleaving) and you case was non consulting the original data layout. Those are now fixed and new versions released so please update both GMT and RemoteS

However, as I guessed above, your image is mostly all blacks. I suggest that you follow exactly that [tutorial](https://www.generic-mapping-tools.org/GMTjl_doc/tutorials/Landsat8/histogram_stretch/) I mentioned above. Just replace the file names for your owns.

Your image

 ![GMTjl_tmp](https://global.discourse-cdn.com/julialang/original/3X/3/a/3a5c748379452440b6ecd800be451fc2393f664d.png)

---

<div class="post-metadata">

**Author:** ![lazarusA](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lazarusa/32/6571_2.png) [@lazarusA](https://discourse.julialang.org/u/lazarusA)\
**Post date:** [December 29, 2022, 10:24pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/13 "2022-12-29T22:24:15Z")

</div>

could you please point me out to the source code for the function that applies the stretching? It looks to do a better job than the ones defined [here](https://juliaimages.org/stable/examples/spatial_transformation/histogram_equalization/#Histogram-equalisation) Any chance to have it in pure images package?

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [December 29, 2022, 11:39pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/14 "2022-12-29T23:39:58Z")

</div>

Sure, it’s in the GMT’s [histogram](https://www.generic-mapping-tools.org/GMTjl_doc/documentation/modules/histogram/index.html#histogram) module. Functions starting [here](https://github.com/GenericMappingTools/GMT.jl/blob/master/src/pshistogram.jl#L274)

In a quick look it doesn’t seems to depend on other GMT code so it should be relatively easy to extract that functionality. But note that it only tries to do a clever pick of the histogram limits to apply the stretch and no histogram equalization. It’s a function specially tailored for UInt16 satellite data. The tutorial I pointed above has a nice example of what it does.

---

<div class="post-metadata">

**Author:** ![lazarusA](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lazarusa/32/6571_2.png) [@lazarusA](https://discourse.julialang.org/u/lazarusA)\
**Post date:** [December 30, 2022, 9:25am UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/15 "2022-12-30T09:25:13Z")

</div>

Yes, I’m interested just in the stretching part. I also have a lot of satellite imagery, unfortunately I had always have problems with the GMT package in my computer. But, this function in specific looks to do actually a really good job. Hence the question.

Edit: They also have [LinearStretching](https://github.com/JuliaImages/ImageContrastAdjustment.jl/pull/28) maybe that one will do the trick for me.

---

<div class="post-metadata">

**Author:** ![EmanuelCastanho](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/emanuelcastanho/32/45245_2.png) [@EmanuelCastanho](https://discourse.julialang.org/u/EmanuelCastanho)\
**Post date:** [December 30, 2022, 11:17am UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/16 "2022-12-30T11:17:11Z")

</div>

I updated `RemoteS` and `GMT` to versions `v0.2.16` and `v0.44.2`, respectively.

```julia
using RemoteS, GMT

Ir = gmtread("Red.tif");
Ir_img = imagesc(Ir);
Ig = gmtread("Green.tif");
Ig_img = imagesc(Ig);
Ib = gmtread("Blue.tif");
Ib_img = imagesc(Ib);

Irgb = truecolor(Ir_img, Ig_img, Ib_img);
imshow(Irgb)

```

This was one solution, I still need to check others suggestions. Thank you!

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [December 30, 2022, 12:22pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/17 "2022-12-30T12:22:21Z")

</div>

> [@lazarusA](#):
>
> unfortunately I had always have problems with the GMT package in my computer.

I guess you are on Linux. Yes, the `libstdc++` compatibility breaking has plagued the automatic installation for several months but it’s now solved. Also, recently `libnetcdf` and `libcurl` started to conflict in recent Linux distros. There is an [open issue](https://github.com/GenericMappingTools/GMT.jl/issues/1052) about it. If it is something else then please open issues.

> [@lazarusA](#):
>
> I also have a lot of satellite imagery

I’m sure the [RemoteS](https://github.com/GenericMappingTools/RemoteS.jl) package would gain a lot from your expertise.

---

<div class="post-metadata">

**Author:** ![EmanuelCastanho](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/emanuelcastanho/32/45245_2.png) [@EmanuelCastanho](https://discourse.julialang.org/u/EmanuelCastanho)\
**Post date:** [December 30, 2022, 4:47pm UTC](https://discourse.julialang.org/t/how-to-create-a-rgb-plot-from-three-different-tif-files/92001/18 "2022-12-30T16:47:47Z")

</div>

Your method also works for my example after some adjustments. However, when I try to apply some transparency to NaN pixels I get a black boundary line. Do you know how this can be removed without changing to white?

Here is the code if someone wants to use in the future:

```julia
# Packages
import ArchGDAL
import Images

# Read example bands
RedRaster = ArchGDAL.read("Red.tif");
RedRasterBand = ArchGDAL.getband(RedRaster, 1);
RedRasterData = Matrix(RedRasterBand);

GreenRaster = ArchGDAL.read("Green.tif");
GreenRasterBand = ArchGDAL.getband(GreenRaster, 1);
GreenRasterData = Matrix(GreenRasterBand);

BlueRaster = ArchGDAL.read("Blue.tif");
BlueRasterBand = ArchGDAL.getband(BlueRaster, 1);
BlueRasterData = Matrix(BlueRasterBand); 

# Construct RGBA
RGBAimg = Matrix(Images.colorview(Images.RGBA, RedRasterData', GreenRasterData', BlueRasterData'));

# Apply transparency to NaNs - I think this can be improved
for i in 1:size(RGBAimg[:])[1]
    if isequal(RGBAimg[i], Images.RGBA{Float32}(NaN32,NaN32,NaN32,1.0f0))
        TempRGBAimg = RGBAimg[i]
        RGBAimg[i] = Images.RGBA{Float32}(TempRGBAimg.r, TempRGBAimg.g, TempRGBAimg.b, 0)
    else # Just to highlight black border
        TempRGBAimg = RGBAimg[i]
        RGBAimg[i] = Images.RGBA{Float32}(1, 0, 0, 1)   
    end
end

# Resize image if needed
SizeF = 10
RGBAimgr = Images.imresize(RGBAimg, SizeF.*size(RGBAimg))

```

 ![Screenshot 2022-12-30 at 15.46.04](https://global.discourse-cdn.com/julialang/original/3X/8/2/8267abbff739e5ab4ec037663f3b20ca34d208dd.png)