wisescootering / infrareddrone

Aerial infrared photography
GNU Affero General Public License v3.0
15 stars 4 forks source link

Make ODM work in multispectral mode #25

Closed balthazarneveu closed 2 years ago

balthazarneveu commented 2 years ago

There is an expected support for multispectral images in ODM as per ODM's doc and as from PR-1057 To enable the mode, it seems like XMP tags are required, A good way would be to use exif tool + custom configuration ( see for instance using config files in exiftool

Specific XMP tags need to be in the images exif Camera:BandName.

Please see this thread multispectral photo

=========================== Update: After testing with the micasense samples provided in the odmdata, I was able to run ODM throug this samples with enabled multi-spectral process and the orthophoto has the 6 bands showing in QGIS.

TODO

To know if you prepared your data correctly you should try to spot this line in the ODM logs [INFO] Reconstruction will use 42 images from auto band

Option : skip_band_alignment may allow to disable another pointless homography estimation in our case where everything has already been aligned run_opensfm.py line 167

<Camera:BandName>Blue</Camera:BandName>
<Camera:CentralWavelength>475</Camera:CentralWavelength>
balthazarneveu commented 2 years ago

Special Exiftool config to write the dedicated XMP tags!

#------------------------------------------------------------------------------
# File:         pix4d.config
#
# Description:  This config file contains tag definitions needed to be able
#               to write Pix4D XMP-Camera tags
#
# Usage:        exiftool -config pix4d.config -XMP-camera:TAG=VAL ...
#
# Requires:     ExifTool version 7.00 or later
#
# References:   1) https://support.pix4d.com/hc/en-us/articles/360016450032-Specifications-of-xmpcamera-tags
#
# Revisions:    2017/12/08 - P. Harvey Created
#               2020/04/02 - PH Updated to current specification
#------------------------------------------------------------------------------

%Image::ExifTool::UserDefined = (
    'Image::ExifTool::XMP::Main' => {
        Camera => { 
            SubDirectory => {
                TagTable => 'Image::ExifTool::UserDefined::Camera',
            },
        },
    },
); 

%Image::ExifTool::UserDefined::Camera = (
    GROUPS => { 0 => 'XMP', 1 => 'XMP-Camera', 2 => 'Camera' },
    NAMESPACE => { 'Camera' => 'http://pix4d.com/camera/1.0/'  },
    WRITABLE => 'string',
    Yaw             => { Writable => 'real' },
    Pitch           => { Writable => 'real' },
    Roll            => { Writable => 'real' },
    IMUSampleSize   => { Writable => 'integer' },
    IMUTimeOffset   => { Writable => 'integer' },
    LineReadoutTime => { Writable => 'integer' },
    IMUFrequency    => { Writable => 'real' },
    PrincipalPoint  => { },
    ModelType       => { },
    PerspectiveFocalLength => { Writable => 'real' },
    PerspectiveDistortion  => { },
    IMULinearVelocity => { },
    GPSXYAccuracy   => { Writable => 'real' },
    GPSZAccuracy    => { Writable => 'real' },
    FlightUUID      => { },
    CentralWaveLength => { },
    BandName        => { List => 'Seq' },
    RigName         => { },
    RigCameraIndex  => { Writable => 'integer' },
    BandName        => { List => 'Seq' },
    IMUAngularVelocity => {
        Binary => 1,
        ValueConv => 'Image::ExifTool::XMP::DecodeBase64($val)',
        ValueConvInv => 'Image::ExifTool::XMP::EncodeBase64($val)',
    },
    # added 2020/04/02 (ref 1)
    FisheyeAffineMatrix     => { },
    FisheyeAffineSymmetric  => { },
    FisheyePolynomial       => { },
    RigRelatives            => { },
    PerspectiveFocalLengthUnits => { },
    CaptureUUID             => { },
    CentralWavelength       => { List => 'Seq' },
    WavelengthFWHM          => { List => 'Seq' },
    BlackCurrent            => { List => 'Seq' },
    BandSensitivity         => { List => 'Seq' },
    SunSensor               => { List => 'Seq' },
    SunSensorExposureTime   => { List => 'Seq' },
    SunSensorSensitivity    => { List => 'Seq' },
    InvalidPixel            => { List => 'Seq' },
    VignettingPolynomial    => { List => 'Seq' },
    VignettingCenter        => { List => 'Seq' },
    VignettingPolynomial2DName => { List => 'Seq' },
    VignettingPolynomial2D  => { List => 'Seq' },
    ColorTransform          => { List => 'Seq' },
    IsNormalized            => { },
    Albedo                  => { List => 'Seq' },
    ReflectArea             => { List => 'Seq' },
    CalibrationPicture      => { Writable => 'integer' },
    GyroRate                => { Writable => 'real' },
    IMUPitchAccuracy        => { Writable => 'real' },
    IMURollAccuracy         => { Writable => 'real' },
    IMUYawAccuracy          => { Writable => 'real' },
    NominalCameraDistance   => { Writable => 'real' },
    AboveGroundAltitude     => { Writable => 'real' },
    SunSensorYaw            => { Writable => 'real' },
    SunSensorPitch          => { Writable => 'real' },
    SunSensorRoll           => { Writable => 'real' },
    SunSensorRelativeRotation => { Writable => 'real', List => 'Seq' },
    TransformAlpha          => { List => 'Seq' },
    TransformBeta           => { List => 'Seq' },
    TransformGamma          => { List => 'Seq' },
    SensorBitDepth          => { Writable => 'integer' },
    SensorTemperature       => { Writable => 'real' },
    # (ref https://community.pix4d.com/t/camera-sun-irradiance-and-sun-angle-in-red-text/3290)
    IrradianceRelativeRotation => { },
);

1;  #end
balthazarneveu commented 2 years ago

Sample

<?xpacket begin="´╗┐" id="W5M0MpCehiHzreSzNTczkc9d"?>
<x:xmpmeta xmlns:x="adobe:ns:meta/" x:xmptk="XMP Core 4.4.0">
   <rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#">
      <rdf:Description rdf:about="Pix4D Camera Information"
            xmlns:Camera="http://pix4d.com/camera/1.0">
         <Camera:RigName>Altum</Camera:RigName>
         <Camera:BandName>Blue</Camera:BandName>
         <Camera:CentralWavelength>475</Camera:CentralWavelength>
         <Camera:WavelengthFWHM>20</Camera:WavelengthFWHM>
         <Camera:ModelType>perspective</Camera:ModelType>
         <Camera:PrincipalPoint>3.57634,2.47155</Camera:PrincipalPoint>
         <Camera:PerspectiveFocalLength>7.7301952350000001</Camera:PerspectiveFocalLength>
         <Camera:PerspectiveFocalLengthUnits>mm</Camera:PerspectiveFocalLengthUnits>      
         <Camera:PerspectiveDistortion>
            <rdf:Seq>
               <rdf:li>-0.13225690000000001</rdf:li>
               <rdf:li>0.2220538</rdf:li>
               <rdf:li>-0.13513149999999999</rdf:li>
               <rdf:li>0.00010770529999999999</rdf:li>
               <rdf:li>1.3645050000000001e-05</rdf:li>
            </rdf:Seq>
         </Camera:PerspectiveDistortion>
         <Camera:VignettingCenter>
            <rdf:Seq>
               <rdf:li>1000.2918</rdf:li>
               <rdf:li>770.69569999999999</rdf:li>
            </rdf:Seq>
         </Camera:VignettingCenter>
         <Camera:VignettingPolynomial>
            <rdf:Seq>
               <rdf:li>0.00012479480000000001</rdf:li>
               <rdf:li>-1.195792e-06</rdf:li>
               <rdf:li>3.0344440000000002e-09</rdf:li>
               <rdf:li>-5.4909070000000004e-12</rdf:li>
               <rdf:li>4.3365830000000002e-15</rdf:li>
               <rdf:li>-1.215311e-18</rdf:li>
            </rdf:Seq>
         </Camera:VignettingPolynomial>
         <Camera:BandSensitivity>0.075946808165865437</Camera:BandSensitivity>
         <Camera:RigRelatives>0.098132, 0.128936, -0.162253</Camera:RigRelatives>
         <Camera:RigCameraIndex>0</Camera:RigCameraIndex>
         <Camera:Irradiance>96.617247953162277</Camera:Irradiance>
         <Camera:IrradianceYaw>-106.26867215486583</Camera:IrradianceYaw>
         <Camera:IrradiancePitch>-10.510424606891545</Camera:IrradiancePitch>
         <Camera:IrradianceRoll>1.0343961122063423</Camera:IrradianceRoll>
      </rdf:Description>
      <rdf:Description rdf:about="Pix4D Camera Information"
            xmlns:MicaSense="http://micasense.com/MicaSense/1.0/">
         <MicaSense:BootTimestamp>126</MicaSense:BootTimestamp>
         <MicaSense:RadiometricCalibration>
            <rdf:Seq>
               <rdf:li>0.00020720760000000001</rdf:li>
               <rdf:li>-1.516716e-10</rdf:li>
               <rdf:li>3.7712690000000002e-05</rdf:li>
            </rdf:Seq>
         </MicaSense:RadiometricCalibration>
         <MicaSense:ImagerTemperatureC>26.219999999999999</MicaSense:ImagerTemperatureC>
         <MicaSense:FlightId>XUhEh4aLKBr7GytxuWs5</MicaSense:FlightId>
         <MicaSense:CaptureId>wNmAilHDCTQ1BacA1EAF</MicaSense:CaptureId>
         <MicaSense:TriggerMethod>0</MicaSense:TriggerMethod>
         <MicaSense:PressureAlt>0</MicaSense:PressureAlt>
         <MicaSense:DarkRowValue>
            <rdf:Seq>
               <rdf:li>3832</rdf:li>
               <rdf:li>3832</rdf:li>
               <rdf:li>3832</rdf:li>
               <rdf:li>3832</rdf:li>
            </rdf:Seq>
         </MicaSense:DarkRowValue>
      </rdf:Description>
      <rdf:Description rdf:about="Pix4D Camera Information"
            xmlns:DLS="http://micasense.com/DLS/1.0/">
         <DLS:Serial>None</DLS:Serial>
         <DLS:SwVersion>None</DLS:SwVersion>
         <DLS:CenterWavelength>475</DLS:CenterWavelength>
         <DLS:Bandwidth>20</DLS:Bandwidth>
         <DLS:TimeStamp>138</DLS:TimeStamp>
         <DLS:SpectralIrradiance>96.617247953162277</DLS:SpectralIrradiance>
         <DLS:HorizontalIrradiance>-118.25894641637043</DLS:HorizontalIrradiance>
         <DLS:DirectIrradiance>146.57127778962638</DLS:DirectIrradiance>
         <DLS:ScatteredIrradiance>7.0271605307142071</DLS:ScatteredIrradiance>
         <DLS:SolarElevation>4.1667182584896265</DLS:SolarElevation>
         <DLS:SolarAzimuth>0.075290360150873717</DLS:SolarAzimuth>
         <DLS:EstimatedDirectLightVector>
            <rdf:Seq>
               <rdf:li>0.27044735266290615</rdf:li>
               <rdf:li>-0.69450682803864394</rdf:li>
               <rdf:li>-0.66672220245416103</rdf:li>
            </rdf:Seq>
         </DLS:EstimatedDirectLightVector>
         <DLS:Yaw>-1.854738220824826</DLS:Yaw>
         <DLS:Pitch>-0.18344151517288818</DLS:Pitch>
         <DLS:Roll>0.01805361792782938</DLS:Roll>
      </rdf:Description>
   </rdf:RDF>
</x:xmpmeta>
balthazarneveu commented 2 years ago

When processing the sample, you get

[INFO]    Finished merge stage
[INFO]    Running opensfm stage
[INFO]    Reconstruction will use 42 images from auto band

During the final stages, you'll get extra logs regarding multispectral image alignment such as

 Can't use features matching, will use ECC (this might take a bit)
[INFO]    Pyramid levels: 2
[INFO]    Computing ECC pyramid level 0
[INFO]    IMG_0188_4.tif --> IMG_0188_1.tif good match
[INFO]    Computing ECC pyramid level 1

image

balthazarneveu commented 2 years ago

image Here's the Micasense samples processed by ODM then displayed in QGIS

image

balthazarneveu commented 2 years ago

Sample code to write 4 single channels TIF to mimick multispectral

import irdrone.process as pr
from pathlib import Path
import numpy as np

dir = Path(r"FLY-test")
dir_vis = dir/"mapping_VIS"/"images"
dir_nir = dir/"mapping_NIR_local"/"images"
out_dir = dir/ "mapping_MS" / "images"
out_dir.mkdir(parents=True, exist_ok=True)

for (idx, img_pth) in enumerate(sorted(list(dir_vis.glob("*.jpg")))):
    img = pr.Image(img_pth)
    img.get_data()
    img_nir_pth = dir_nir/img_pth.name.replace("VIS", "NIR_local_")
    assert img_nir_pth.exists()
    nir_img = pr.Image(img_nir_pth)
    nir_img.get_data()

    ms_data = np.zeros((img._data.shape[0], img._data.shape[1], 4))
    ms_data[:, :, :3] = img._data
    ms_data[:, :, 3] = nir_img._data[:, :, 1]
    img._data = ms_data
    ms_data.shape
    out_name = f"{(idx+1):04d}"
    img.save_multispectral(out_dir/out_name)

image

After adding custom XMP support, XNview detects the tag correctly & the XMP tags are similar to the ones in a Micasense TIF image

Visible RGB band in QGIS image NIR band in QGIS - perfectly aligned with visible image

balthazarneveu commented 2 years ago

image VS image

balthazarneveu commented 2 years ago

image 👍 WebODM supports our multispectral TIFs 👍

balthazarneveu commented 2 years ago

QGIS_multispectral_maps_creation.pdf Documentation by @AlainNeveu on how to create VIR maps & NDVI maps directly in QGIS

balthazarneveu commented 2 years ago

@AlainNeveu @Neveu-Alain
Please put the

AlainNeveu commented 2 years ago

Echelle de couleur NDVI:

Créer un fichier ...Air-Mission\NDVI-00.txt (par exemple dans le dossier Air-Mission

# Fichier d'export QGIS de palette de couleur
INTERPOLATION:INTERPOLATED
-0.94999999999999996,0,0,0,255,-0,9500
-0.83009188361408881,1,1,1,255,-0,8301
-0.70699999999999996,140,140,138,255,-0,7070
-0.59299999999999997,180,141,92,255,-0,5930
-0.49473684210526309,190,185,88,255,-0,4947
-0.41684782608695647,166,214,102,255,-0,4168
-0.34399999999999986,122,246,34,255,-0,3440
-0.29299999999999993,102,235,7,255,-0,2930
-0.22999999999999987,11,184,8,255,-0,2300
-0.16699999999999993,8,140,52,255,-0,1670
-0.06956521739130428,2,102,30,255,-0,0696
0.05835889570552166,0,91,92,255,0,0584
0.2483695652173914,62,0,88,255,0,2484
0.40000000000000013,126,0,130,255,0,4000

Dans le projet QGIS sélectionner la couche NDVI puis dbl clic > onglet Symbologie > type de rendu "pseudo-couleur à bande unique" Clic sur icone en forme de dossier "charger une palette de couleur depuis un fichier" . Sélectionner ... Air-Mission\NDVI-00.txt

Il est possible d'éditer une palette et de la modifier très simplement: palette de couleur > éditer la palette de couleurs ... On peut aussi la sauvegarder pour la faire apparaître dans la liste des palettes couleur de QGIS palette de couleur > sauvegarder la palette de couleurs > Nom : NDVI (cocher ajouter aux marque page)

NDVI-00.txt