AnalyzeFieldSpectrometerMeasurements/ApplyLUTs.m at master · JorisTimmermans/AnalyzeFieldSpectrometerMeasurements · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
clear all
close all
fclose all;
clc

%% define datadir
dirname         = SelectDirectory;

%% Read Canopy spectrum
files           =   dir([dirname,'/*.sed']);
for j=1:length(files)
    filename    =   files(j).name;
    fid         =   fopen([dirname,'/',filename],'r');
    data        =   textscan(fid,'%f %f','headerlines',27);
    fclose(fid);

    wl          =   data{1};
    Refl(:,j)   =   data{2}/100;
end

%filter out atmospheric effects
iatm            =   (wl>1800 & wl<2000) | (wl>2300);
Refl(iatm,:)    =   NaN;

%% load LUT data
lut_dir         =   [dirname,'\LUT\'];
if ~exist([lut_dir,'LUT_SLC.mat'],'file')
    warndlg('No LUT file found, please run "CreateLUTs.m" ')
    return
end
load([lut_dir,'LUT_SLC.mat'],'BandsReflectance','Vegpar','Leafbrown','Leafgreen')

%% plot data
subplot(2,1,1)
plot(wl,Refl)
title('Measurements')
subplot(2,1,2)
plot(wl,BandsReflectance)
title('LUT')

nLUT = size(BandsReflectance,2);
nwl             =   length(wl);
filter          =   ones(nwl,1);
for j=1:size(Refl,2)
    can_refl    =   Refl(:,j);
    Can_Refl    =   can_refl(:, ones(1,nLUT))*pi;


%     plot(wl,BandsReflectance,'b',wl,can_refl,'r')
%     keyboard
    E           =   BandsReflectance - Can_Refl;

    RMSE        =   sqrt(nanmean(E(:,:).^2));
    [~,imin]    =   min(RMSE);

%     rMSE        =   sqrt(nanmean(E./Can_Refl.^2));
%     [~,imin]    =   min(RMSE);

    for ilut=1:size(BandsReflectance,2)
        c           =   corrcoef(BandsReflectance(50:1440,ilut), can_refl(50:1440,1));
        C(ilut)     =   c(1,2);
    end
    [~,imin]    =   max(C);

%     keyboard
    leafbrown   =   Leafbrown(:,imin);
    leafgreen   =   Leafgreen(:,imin);
    vegpar      =   Vegpar(:,imin);

    lai(j)         =   vegpar(1);
%     lidfa       =   vegpar(2);
%     lidfb       =   vegpar(3);
%     hot         =   vegpar(4);
%     fbrown      =   vegpar(5);
%     diss        =   vegpar(6);
%     clumping    =   vegpar(7);
%     zeta        =   vegpar(8);

    cab(j)         =   leafgreen(1);
    cw(j)          =   leafgreen(2);
    cdm(j)         =   leafgreen(3);
    cs(j)          =   leafgreen(4);
    N(j)           =   leafgreen(5);
end

%%
figure
subplot(3,1,1)
plot(cab)
ylabel('cab')
subplot(3,1,2)
plot(lai)
ylabel('lai')
subplot(3,1,3)
plot(cw)
ylabel('cw')
% save([lut_dir,'LUT_SLC.mat'],'Leafgreen','Leafbrown','Vegpar','BandsReflectance','wl')

% save(