Converting pixel intensity to temperature from TIFF image - Grayscale values much lower than expected

Hello all,

I am working on a matlab script to input a TIFF file that was converted from RAW using the open source program dcraw. I am then loading the TIFF into the matlab script and using some conversions, getting temperature from intensity. There are three methods in the code Green/Red intensity ratio, grayscale intensity, and a hybrid that combines both. My issue lies in the grayscale values, which are derived from the RGB values, being much lower than the G/R ratio results when they should be relatively close. I am not sure what is wrong or why this is happening. Perhaps I am not handling the grayscale math correctly or if there is something wrong with the RAW->TIFF conversion.

The script is here:

NOTE: unfortunately the pixel->temperature conversion equations are proprietary and I cannot share them but they are from peer-reviewed literature so they are likely not the cause of the discrepancy. If no other issues are found I may look into those equations specifically though.

clear all
tic
format long g
myFolder = 'Your directory with the TIFF';
% Get a list of all files in the folder with the desired file name pattern.
filePattern = fullfile(myFolder, '*.tiff'); % Change to whatever pattern you need.
theFiles = dir(filePattern);

for k = 1: length(theFiles)
    photoname = theFiles(k).name;
    fullFileName = fullfile(theFiles(k).folder, photoname);
    
    A = imread(fullFileName);
    rdc = 7.7;
    gdc = 11.5;
    bdc = 4.4;
    GSdc = (rdc+gdc+bdc)/3;
    iso = 6400;
    t = 0.2;
    f = 2.4;
    y1 = 2000; % Entire picture is 5496 by 3672 pixels.
    y2 = 2150;
    x1 = 1800;
    x2 = 2200;
    nrows = y2 - y1;
    ncolumns = x2 - x1;
    rows = linspace(y1,y1+nrows-1,nrows);
    for i = 1:ncolumns
        c(i,1:nrows) = linspace(x1+i-1,x1+i-1,nrows);
    end

    for i = 1:ncolumns
        RGB(1:nrows,(((4*i)-3):((4*i)-1))) = impixel(A,c(i,1:nrows),rows);
    end
    %% Separating RGB columns
    rawred = RGB(1:end,1:4:end);
    rawred2 = zeros(nrows,ncolumns);
    for i = 1:nrows
        for j = 1:ncolumns
            if rawred(i,j) > 65534
                rawred2(i,j) = 0;
            else
                rawred2(i,j) = rawred(i,j);
            end
        end
    end
    rawgreen = RGB(1:end,2:4:end);
    rawblue = RGB(1:end,3:4:end);
    %% Normalized Intensities
    rawGS = (rawred2+rawgreen+rawblue)/3;
    normred = ((rawred2-rdc)*(f^2))/(iso*t);
    normgreen = ((rawgreen-gdc)*(f^2))/(iso*t);
    normblue = ((rawblue-bdc)*(f^2))/(iso*t);
    normGS = ((rawGS-GSdc)*(f^2))/(iso*t);
    %% Finding log of pixel G/R ratio to use for ratio pyrometry curve fit
    GR = normgreen./normred;
    GR2 = log10(GR);
    IP = imag(GR2);
    logGR = zeros(nrows,ncolumns);
    for i=1:nrows
        for j=1:ncolumns
            if IP(i,j) > 0
                logGR(i,j) = 0;
            else
                logGR(i,j) = GR2(i,j);
            end
        end
    end
    %% Finding log of pixel Grayscale ratio 
    GS2 = log10(normGS);
    IPGS = imag(GS2);
    logGS = zeros(nrows,ncolumns);
    for i=1:nrows
        for j=1:ncolumns
            if IPGS(i,j) > 0
                logGS(i,j) = 0;
            else
                logGS(i,j) = GS2(i,j);
            end
        end
    end
    %THIS SECTION HAS BEEN CHANGED DUE TO RESTRICTIONS ON RELEASE OF
    %PROPRIETARY INFORMATION
    %% Application of ratio pyrometry curve fit
    RT = (1*(logGR.^3) + 1*(logGR.^2) + 1*(logGR) + 1); % New Camera Ratio Temp
    %% Grayscale Pyrometry Fit
    GST = (1.*(logGS.^2))+(1.*logGS)+1;
    %% Hybrid Pyrometry Method
    HT = 1
    %% Removal of 600-1200 range (Ratio)
    %END OF SECTION
    RT2 = zeros(nrows,ncolumns);
    for i=1:nrows
        for j=1:ncolumns
            if RT(i,j) < 600
                RT2(i,j) = 0;
            elseif RT(i,j) > 1200
                RT2(i,j) = 0;
            else
                RT2(i,j) = RT(i,j);
            end
        end
    end
    %% Removal of 600-1200 range (Grayscale)
    GST2 = zeros(nrows,ncolumns);
    for i=1:nrows
        for j=1:ncolumns
            if GST(i,j) < 600
                GST2(i,j) = 0;
            elseif GST(i,j) > 1200
                GST2(i,j) = 0;
            else
                GST2(i,j) = GST(i,j);
            end
        end
    end
    %% Removal of 600-1200 range (Hybrid)
    HT2 = zeros(nrows,ncolumns);
    for i=1:nrows
        for j=1:ncolumns
            if HT(i,j) < 600
                HT2(i,j) = 0;
            elseif HT(i,j) > 1200
                HT2(i,j) = 0;
            else
                HT2(i,j) = HT(i,j);
            end
        end
    end
    %% 50% rule (Ratio)
    RT3 = zeros(nrows,ncolumns);
    for i = 4:(nrows-4)
        for j = 4:(ncolumns-4)
            if nnz(RT2(i-3:i+3,j-3:j+3)) < 25
                RT3(i,j) = 0;
            else
                RT3(i,j) = RT2(i,j);
            end
        end
    end
    %% 50% rule (Grayscale)
    GST3 = zeros(nrows,ncolumns);
    for i = 4:(nrows-4)
        for j = 4:(ncolumns-4)
            if nnz(GST2(i-3:i+3,j-3:j+3)) < 25
                GST3(i,j) = 0;
            else
                GST3(i,j) = GST2(i,j);
            end
        end
    end
    %% 50% rule (Hybrid)
    HT3 = zeros(nrows,ncolumns);
    for i = 4:(nrows-4)
        for j = 4:(ncolumns-4)
            if nnz(HT2(i-3:i+3,j-3:j+3)) < 25
                HT3(i,j) = 0;
            else
                HT3(i,j) = HT2(i,j);
            end
        end
    end
    %% MATLAB display
    %RT4 = flip(RT3);
    RT4=RT3;
    GST4=GST3;
    HT4=HT3;
    RATraw = mean(nonzeros(RT4));
    RatioAverageTemp = RATraw(~isnan(RATraw))
    GSATraw = mean(nonzeros(GST4));
    GrayscaleAverageTemp = GSATraw(~isnan(GSATraw))
    HATraw = mean(nonzeros(HT4));
    HybridAverageTemp = HATraw(~isnan(HATraw))
        figure(1)
        P = imagesc(RT4,[600 1200]);
        shading interp
        colormap jet
        colorbar;
        title([photoname ' Green-Red Ratio Analysis']);
        xlabel(['Ratio Average Temperature = ' num2str(RatioAverageTemp)]);
            figure(2)
            Q = imagesc(GST4,[600 1200]);
            shading interp
            colormap jet
            colorbar;
            title([photoname ' Grayscale Analysis']);
            xlabel(['Grayscale Average Temperature = ' num2str(GrayscaleAverageTemp)]);
                figure(3)
                R = imagesc(HT4,[600 1200]);
                shading interp
                colormap jet
                colorbar;
                title([photoname ' Hybrid Analysis']);
                xlabel(['Hybrid Average Temperature = ' num2str(HybridAverageTemp)]);
    toc
    %% Statistics
    RSTD = std(nonzeros(RT4));
    GSSTD = std(nonzeros(GST4));
    HSTD = std(nonzeros(HT4));
    RSNR = RatioAverageTemp/RSTD
    GSSNR = GrayscaleAverageTemp/GSSTD
    HSNR = HybridAverageTemp/HSTD
    RatioCI95 = (1.96*RSTD)/sqrt(nnz(RT4))
    GrayscaleCI95 = (1.96*GSSTD)/sqrt(nnz(GST4))
    HybrdiCI95 = (1.96*HSTD)/sqrt(nnz(HT4))
end

Here are the RAW and TIFF images. For those familiar with dcraw the input command used was

dcraw -v -w -H 0 -o 1 -q 3 -4 -T

in command line using dcraw v9.27

RAW IMAGE

TIFF Conversion

Hopefully the google drive links dont change the files from how they are stored originally.

Answer 1

[Image Analyst]

"My issue lies in the grayscale values, which are derived from the RGB values, being much lower than the G/R ratio results when they should be relatively close."

Why do you think grayscale values, presumably a weight sum of R, G, and B like from rgb2gray() function, should be close to the G/R ratio? I mean gray scale values are likely in the range 0-255 while the ratio can be something less than 1 or as much as infinity (if R=0). Typically I'd think the ratio would often be in the 0.5 to 2 or 3 range while gray scales could be 50 or 150 or 200 or so. Why should there be any reason for the values to be close?

[Mahdi Tlemsani]

Hello, Thank you for taking the time to help me with this question.

The object of the script is to take a raw image, convert to a 16-bit TIFF image which has intensity values 0-65535, normalize them, and use an equation derived from testing to convert those intensities to temperature values. The technique is referred to in research papers as color pyrometry I believe. My results are showing when using an equation that take G/R as input my temperatures are appropriate but the equation meant to use grayscale intensities does not produce correct results. The equation has been used succesfully in previous published literature so it is not the issue so I figured it had to be something with my technique of either converting from RAW to TIFF in dcraw or my use of calculating those values in the script shown here.

Hopefully this clears things up.

[Image Analyst]

Well, what can I say? All I can say is to make sure your code follows the formula(s) in your paper.

[DGM]

Were they operating on sRGB or linear image data?

[Walter Roberson]

rgb2gray gives an idea of how bright the pixel is. G/R gives an idea of how pure gray it is.

[Image Analyst]

I don't see how a visible RGB digital camera would give you temperature unless the object was SO hot that it was literally glowing in the visible, like red hot or white hot so you could use Planck's law for black body radiation. You're talking like 2000 to 15000 degrees Kelvin - extremely hot, nowhere near room temperature.

[Walter Roberson]

According to the metadata, a Sony RX10M3 was used, which is one of Sony's compact digital cameras. https://www.sony.ca/en/electronics/cyber-shot-compact-cameras/dsc-rx10m3

The metadata also says that the shutter speed was 1/5 second, at ISO 6400 (F2.4), sRGB IEC61966-2.1.

That's a really long shutter speed for ISO 6400! Equivalent light exposure to two minutes at ISO 100. And yet the images are certainly not saturated.

So it does appear to be a visible light red glow... of something that is not glowing much.

... Ah, I see sites saying that you can take IR photos with the RX10M3 and RX10M4 if you use a Hoya RM72 Infrared Filter and increase the exposure... so maybe it is NIR.

[Mahdi Tlemsani]

To address some of the comments.

I believe the conversion in dcraw used sRGB for this image. I should have specified sooner but the image is of a glowing red hot ember as it passes by the camera at 6-8m/s. We are looking at a relatively novel technique called color camera pyrometry instead of using IR sensors.

For the pixel value->temp formulas im relatively sure they are ok. The issue is my grayscale intensity values seem low in the image and I could not figure out why.

[Walter Roberson]

@Mahdi Tlemsani

The .ARW image is 3648 x 5472, but the dcraw is 3672 x 5496

If you use Sony's Image Edge program to convert the ARW to TIF, and use imshowpair() to compare against the dcraw (rotated by 270 degrees), you will see that they do not quite align. If you suppose that one might have padding along an unfortunate edge, and so flip the images upside down before doing the imshowpair()... then it still does not align. The disagreement varies along the length.

Please confirm that they are in fact intended to be the same image.

[Mahdi Tlemsani]

The raw and tiff images are intended to be exactly the same. I am not sure what is causing this loss of pixels in dcraw other than possibly the interpolation method used. dcraw has very poor documentation but that is the only input I had that might do what you're describing.

[Walter Roberson]

dcraw is adding 24 pixels of height and width.

I wonder if it is padding on all sides?

Let me see if I can find some image registration...

[Walter Roberson]

When I use Sony's Image Edge to export the file to 16 bit TIFF, and use the Image Registration Estimater against the dcraw image, the closest I can get is a 90% match (using Phase Correlation), and the alignment is visibly off, especially towards the ends.

I suggest experimenting with converting the ARW file using Sony's Image Edge; I think it will give you different results.

[Image Analyst]

Can't you just crop off the padding if it's a known size?

If not, use normxcorr2() to find the smaller image in the bigger image. Demo attached.

[Walter Roberson]

Image Registration Editor is saying that it is not just a matter of cropping off padding, that the data in the cdraw file is only about 88% similar to the data in the ARW file.

[Mahdi Tlemsani]

Could this be interpolation used in the conversion? Additionally in my script I also remove some data if the temperature is below 600C pr above 1200C from the temperature equations and I remove some imaginary numbers that occur during the process.

I've also attached the paper I used as a reference for this code so possibly it could help answer the question. In the end my original issue is that grayscale pixel intensities seems off. Perhaps I am simply not calculating them right, (R + G + B)/3

[Walter Roberson]

Gray intensity is not the mean of R, G, B! They have to be weighted. G has the highest weight, then R, and B has a fairly small weight. rgb2gray does this.

The eye is most sensitive to green brightness, least sensitive to blue brightness. (It is, however, most sensitive to blue contrast. Differences in blue are perceived more readily than differences in the others)