diff --git a/.gitattributes b/.gitattributes
new file mode 100644
index 0000000000000000000000000000000000000000..44a043b3aa6964811a2d76e7d3ce982a0a98629f
--- /dev/null
+++ b/.gitattributes
@@ -0,0 +1,2 @@
+*.nc filter=lfs diff=lfs merge=lfs -text
+*.mat filter=lfs diff=lfs merge=lfs -text
diff --git a/RU/.gitattributes b/RU/.gitattributes
new file mode 100644
index 0000000000000000000000000000000000000000..9abd2058819fae6e5f68d95ef7ffd1a2abad7650
--- /dev/null
+++ b/RU/.gitattributes
@@ -0,0 +1 @@
+*.nc filter=lfs diff=lfs merge=lfs -text
diff --git a/RU/Ep_approx.m b/RU/Ep_approx.m
new file mode 100644
index 0000000000000000000000000000000000000000..c8d3d8bb74e19707e86f79eba313c8bfcf8aad2c
--- /dev/null
+++ b/RU/Ep_approx.m
@@ -0,0 +1,52 @@
+function E_p = Ep_approx(prpt,alph,delta_ext,L_cbb,L_hbb,B_ssm,L,d_cbb,d_hbb)
+%
+% function E_p = Ep_approx(prpt,alph,delta_ext,L_cbb,L_hbb,B_ssm,L,d_cbb,d_hbb)
+%
+% The SSM and sensor act like a polarizer-analyzer pair. This function calculates
+% the approximate radiance error generated by the polarizer-analyzer pair
+%
+% Inputs
+% prpt: p_r * p_t for sensor + SSM (unitless) [nchan x 1]
+% alph: analyzer angle alpha for sensor (radians) [scalar]
+% delta_ext: polarizer angle for N scene views (radians) [1 x N]
+% L_cbb: Radiance for cold cal reference (Rad units)) [nChan x N]
+% L_hbb: Radiance for hot cal reference (Rad units)) [nChan x 1]
+% B_ssm: Planck radiance for scene mirror temp (Rad units)) [nChan x 1]
+% L: Scene radiance for N scene views (Rad units)) [nChan x 1] or [nChan x N]
+% d_cbb: polarizer angle for CBB view (radians) [scalar]
+% d_hbb: polarizer angle for HBB view (radians) [scalar]
+%
+% Outputs
+% E_p calibration bias due to polarization (Rad units) [nChan x N]
+%
+% 2020-08-19 JKT
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+
+if size(L,2) == 1
+ L = repmat(L,1,length(delta_ext));
+elseif size(L,2) ~= length(delta_ext)
+ fprintf(1, 'Ep_approx: ERROR Invalid scene radiance dimensions\n');
+end
+if abs(d_cbb) > 2*pi
+ fprintf(1, 'Ep_approx: ERROR delta CBB out of range, should be in units of radians\n')
+ return
+end
+if abs(d_hbb) > 2*pi
+ fprintf('Ep_approx: ERROR delta HBB out of range, should be in units of radians\n')
+ return
+end
+E_p = NaN*ones(length(L_cbb),length(delta_ext)); % N_{d.wn} x N_{d.delta_ext}
+
+for ii = 1:length(delta_ext)
+ d_ext = delta_ext(ii);
+ L_ext = L(:,ii);
+ % approximate equation
+ term1 = prpt.*(L_ext*cos(2*(d_ext - alph)));
+ term2 = -prpt.*(L_hbb.*((L_ext - L_cbb)./(L_hbb - L_cbb))*cos(2*(d_hbb - alph)));
+ term3 = -prpt.*(L_cbb.*((L_hbb - L_ext)./(L_hbb - L_cbb))*cos(2*(d_cbb - alph)));
+ term4 = -B_ssm.*prpt.*(cos(2*(d_ext - alph)));
+ term5 = +B_ssm.*prpt.*(((L_ext - L_cbb)./(L_hbb - L_cbb))*cos(2*(d_hbb - alph)));
+ term6 = +B_ssm.*prpt.*(((L_hbb - L_ext)./(L_hbb - L_cbb))*cos(2*(d_cbb - alph)));
+ E_p(:,ii) = term1 + term2 + term3 + term4 + term5 + term6;
+end
+return
diff --git a/RU/ICTradModelRU_L1b_aux.m b/RU/ICTradModelRU_L1b_aux.m
new file mode 100644
index 0000000000000000000000000000000000000000..9cec1a966b0dcd3c01d5423b1aade72cbffa77cb
--- /dev/null
+++ b/RU/ICTradModelRU_L1b_aux.m
@@ -0,0 +1,262 @@
+function B = ICTradModelRU_L1b_aux(sensor,T_ICT,T_CRIS,ICT_Param,ru_term,ru_value)
+%
+% function B = ICTradModelRU_L1b_aux(sensor,T_ICT,T_CRIS,ICT_Param,ru_term,ru_value)
+%
+% Compute ICT predicted radiance using ICT radiometric model prescribed by ATBD
+% with perturbation to a single term
+%
+% Inputs
+% sensor sparse spectral params structure
+% .band
+% .v
+% .output_range
+% T_ICT ICT Temperature
+% T_CRIS CrIS temperature struct
+% ICT_Param ICT radiometric model parameters
+% ru_term string indicating ru_term to perturb
+% ru_value size of perturbation
+%
+% Outputs
+% B Structure of predicted ICT radiance variables
+% .emitted emitted radiance
+% .wn wavenumber scale
+% .SSM_Baffle SSM Baffle term reflected radiance
+% .ICT_baffle ICT baffle term reflected radiance
+% .OMA OMA term reflected radiance
+% .BS_warm BS_warm term reflected radiance
+% .BS_cold BS_cold term reflected radiance
+% .Space Earth/Space term reflected radiance
+% .reflected total reflected radiance
+% .total total radiance (emitted and reflected)
+% .V view factors (struct)
+% .e emissivities (struct)
+% .r reflectivies (.SSM)
+% .T temperatures (struct)
+% .e_ICT ICT emissivity
+% .T_ICT ICT temperature
+% .ru_term RU term to perturb
+% .ru_value value of perturbation
+% .no_perturb_val
+% .perturbed_val
+%
+% 2020-08-19 JKT
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+
+wn = sensor.v_onaxis;
+
+ICT_wnum = ICT_Param.(sprintf('%s_ICT_wnum',lower(sensor.band)));
+ICT_emiss = ICT_Param.(sprintf('%s_ICT_emiss',lower(sensor.band)));
+ICT_Baffle_emiss = ICT_Param.(sprintf('%s_ICT_Baffle_emiss',lower(sensor.band)));
+Housing_emiss = ICT_Param.(sprintf('%s_Housing_emiss',lower(sensor.band)));
+ScanBaffle_emiss = ICT_Param.(sprintf('%s_ScanBaffle_emiss',lower(sensor.band)));
+ScanMirror_emiss = ICT_Param.(sprintf('%s_ScanMirror_emiss',lower(sensor.band)));
+Earth_emiss = ICT_Param.(sprintf('%s_Earth_emiss',lower(sensor.band)));
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% ICT emissivity interpolation to wavenumber grid
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% Interpolate 4-min-EP emissivity to the on-axis sensor grid.
+%% The wavenumber values corresponding to the emissivity values in the engineering packet are on the normal spectral resolution and on user grid
+e_ICT = interp1(ICT_wnum,ICT_emiss,wn,'nearest','extrap');
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% Temperature, emissivity and view factor init
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%% SSM
+%%----------
+%r.SSM = 1-ICT_Param.Band(iband).Emissivity.ScanMirrorPts;
+r.SSM = 1-ScanMirror_emiss;
+
+%% SSM_Baffle
+%%----------
+%e.SSM_Baffle = ICT_Param.Band(iband).Emissivity.ScanBafflePts;
+e.SSM_Baffle = ScanBaffle_emiss;
+T.SSM_Baffle = T_CRIS.ssm_baffle_temp + T_CRIS.ssm_baffle_corr_temp;
+T.SSM_Baffle(isnan(T.SSM_Baffle)) = 0;
+V.SSM_Baffle = ICT_Param.View_Factor.ScanBaffle;
+
+%% ICT_baffle
+%%----------
+%e.ICT_baffle = ICT_Param.Band(iband).Emissivity.ICT_BafflePts;
+e.ICT_baffle = ICT_Baffle_emiss;
+T.ICT_baffle = T_ICT;
+V.ICT_baffle = ICT_Param.View_Factor.ICTBaffle;
+
+%% OMA + Frame + Other
+%%----------
+%e.OMA = ICT_Param.Band(iband).Emissivity.HousingPts;
+e.OMA = Housing_emiss;
+T.OMA = mean(0.5*T_CRIS.oma_struct_temp_1 + 0.5*T_CRIS.oma_struct_temp_2);
+V.OMA = ICT_Param.View_Factor.OMAandFrame;
+
+%% BS_warm
+%%----------
+e.BS_warm = e.OMA;
+T.BS_warm = T.OMA;
+V.BS_warm = ICT_Param.View_Factor.BeamSplitterWarm;
+
+%% BS_cold
+%%----------
+T.BS_cold = T_ICT;
+V.BS_cold = ICT_Param.View_Factor.BeamSplitterCold;
+
+%% Space/Earth Reflected
+%%----------
+e.Space = Earth_emiss;
+T.Space = ICT_Param.Earth_Temperature;
+V.Space = ICT_Param.View_Factor.Space;
+
+no_perturb_val = []; % initialize, will remain at this value for no perturbation
+perturbed_val = []; % initialize, will remain at this value for no perturbation
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% Perturbation
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+switch ru_term
+ case 'e_ICT'
+ no_perturb_val = e_ICT;
+ e_ICT = e_ICT*(1-ru_value);
+ e_ICT(e_ICT > 1) = 1;
+ e_ICT(e_ICT < 0) = 0;
+ perturbed_val = e_ICT;
+
+ case 'T_ICT'
+ no_perturb_val = T_ICT;
+ T_ICT = T_ICT + ru_value;
+ perturbed_val = T_ICT;
+
+ case 'SSM_scan_mirror_baffle'
+ no_perturb_val = T.SSM_Baffle;
+ T.SSM_Baffle = T.SSM_Baffle + ru_value;
+ perturbed_val = T.SSM_Baffle;
+
+ case 'T_ScanBaffleCorrection'
+ no_perturb_val = T.SSM_Baffle;
+ T.SSM_Baffle = T.SSM_Baffle + ru_value;
+ perturbed_val = T.SSM_Baffle;
+
+ case 'T_ICT_baffle'
+ no_perturb_val = T.ICT_baffle;
+ T.ICT_baffle = T.ICT_baffle + ru_value;
+ perturbed_val = T.ICT_baffle;
+
+ case 'T_OMA' % OMA + Frame + Other
+ no_perturb_val = T.OMA;
+ T.OMA = T.OMA + ru_value;
+ perturbed_val = T.OMA;
+
+ case 'T_Frame' % no longer used
+ no_perturb_val = T.Frame;
+ T.Frame = T.Frame + ru_value;
+ perturbed_val = T.Frame;
+
+ case 'T_BS_warm'
+ no_perturb_val = T.BS_warm;
+ T.BS_warm = T.BS_warm + ru_value;
+ perturbed_val = T.BS_warm;
+
+ case 'T_BS_cold'
+ no_perturb_val = T.BS_cold;
+ T.BS_cold = T.BS_cold + ru_value;
+ perturbed_val = T.BS_cold;
+
+ case 'T_Space'
+ no_perturb_val = T.Space;
+ T.Space = T.Space + ru_value;
+ perturbed_val = T.Space;
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% ICT emission
+% (ref: eq 75 in 474-00032, Rev C)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+B.emitted = e_ICT .* bt2rad(wn,T_ICT);
+B.wn = wn;
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Reflected term 1: SSM Baffle
+% (ref: eq 76 in 474-00032, Rev C)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+if V.SSM_Baffle == 0
+ B.SSM_Baffle = zeros(size(B.wn));
+else
+ B.SSM_Baffle = (1 - e_ICT) .* e.SSM_Baffle * V.SSM_Baffle .* bt2rad(wn,T.SSM_Baffle);
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Reflected term 2: ICT Baffle
+% (ref: eq 78a in 474-00032, Rev C)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+if V.ICT_baffle == 0
+ B.ICT_baffle = zeros(size(B.wn));
+else
+ B.ICT_baffle = (1 - e_ICT) .* e.ICT_baffle * V.ICT_baffle .* bt2rad(wn,T.ICT_baffle);
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Reflected term 3a,b,c: OMA, frame, other, Beamsplitter Warm
+% (ref: eq 77 in 474-00032, Rev C)
+% V = 0.214 for OMA, Other, Frame terms in ATBD (= 0.0590 + 0.0100 + 0.1450)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% OMA
+%%----------
+if V.OMA == 0
+ B.OMA = zeros(size(B.wn));
+else
+ B.OMA = (1 - e_ICT) .* e.OMA * r.SSM * V.OMA .* bt2rad(wn,T.OMA);
+end
+
+%% BS_Warm
+%%----------
+if V.BS_warm == 0
+ B.BS_warm = zeros(size(B.wn));;
+else
+ B.BS_warm = (1 - e_ICT) .* e.BS_warm * r.SSM * V.BS_warm .* bt2rad(wn,T.BS_warm);
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Reflected term 4: Beamsplitter Cold
+% (ref: eq 78b in 474-00032, Rev C)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+if V.BS_cold == 1
+ B.BS_cold = (1 - e_ICT) .* r.SSM^2 .* 0.5 * V.BS_cold .* bt2rad(wn,T.BS_cold); % Consistent with ATBD
+else
+ B.BS_cold = (1 - e_ICT) .* r.SSM^2 .* 0.5 * V.BS_cold .* bt2rad(wn,T.BS_cold); % SNPP
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Reflected term 5: Space/Earth View
+% (ref: eq 78c in 474-00032, Rev C)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+if V.Space == 0
+ B.Space = zeros(size(B.wn));;
+else
+ B.Space = (1 - e_ICT) .* V.Space .* e.Space .* bt2rad(wn,T.Space);
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% total reflected radiance into ICT cavity
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+V.Total = V.SSM_Baffle + V.ICT_baffle + V.OMA + V.BS_warm + V.BS_cold + V.Space; % should total 1, for check only
+B.reflected = B.SSM_Baffle + B.ICT_baffle + B.OMA + B.BS_warm + B.BS_cold + B.Space;
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Total radiance from ICT
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+B.total = B.emitted + B.reflected;
+
+B.wn = wn;
+B.V = V;
+B.e = e;
+B.r = r;
+B.T = T;
+B.e_ICT = e_ICT;
+B.T_ICT = T_ICT;
+B.ru_term = ru_term;
+B.ru_value = ru_value;
+B.no_perturb_val = no_perturb_val;
+B.perturbed_val = perturbed_val;
+
+return
diff --git a/RU/J1_L1B_ru_tool_param_file.nc b/RU/J1_L1B_ru_tool_param_file.nc
new file mode 100644
index 0000000000000000000000000000000000000000..f98e4e87cea5f2d2071d47dc942cc3ad92165dab
--- /dev/null
+++ b/RU/J1_L1B_ru_tool_param_file.nc
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:6c80a8d82ae12d8061fed579c187def89e13b80aa43bce271ea6f261a1ba5528
+size 198316
diff --git a/RU/SNPP_L1B_ru_tool_param_file.nc b/RU/SNPP_L1B_ru_tool_param_file.nc
new file mode 100644
index 0000000000000000000000000000000000000000..e284e0cf6c70e9919fcdcb0e38ab53efe105bc4e
--- /dev/null
+++ b/RU/SNPP_L1B_ru_tool_param_file.nc
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:5a9389260928e98dccc68206c87bde00e12b2db4d1118088ce35084901460884
+size 198316
diff --git a/RU/TEST_INPUT/SNDR.J1.CRIS.20191031T0012.m06.g003.L1B.std.v03_08.G.201204050506.nc b/RU/TEST_INPUT/SNDR.J1.CRIS.20191031T0012.m06.g003.L1B.std.v03_08.G.201204050506.nc
new file mode 100644
index 0000000000000000000000000000000000000000..9d85e52448d8c51d25e7ee0b8aa03da6cb54fae4
--- /dev/null
+++ b/RU/TEST_INPUT/SNDR.J1.CRIS.20191031T0012.m06.g003.L1B.std.v03_08.G.201204050506.nc
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:f008b1e83c3574fe879c5a7bd6d1d7f15c2f05c131a164d00136a7f87434cfa8
+size 175523775
diff --git a/RU/TEST_INPUT/SNDR.SNPP.CRIS.20191031T0012.m06.g003.L1B.std.v03_08.G.201219125854.nc b/RU/TEST_INPUT/SNDR.SNPP.CRIS.20191031T0012.m06.g003.L1B.std.v03_08.G.201219125854.nc
new file mode 100644
index 0000000000000000000000000000000000000000..10e7804c1520a381e435f97f5a1b855b43faff11
--- /dev/null
+++ b/RU/TEST_INPUT/SNDR.SNPP.CRIS.20191031T0012.m06.g003.L1B.std.v03_08.G.201219125854.nc
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:fa6c558c697c15cb52d8bb28305d4c9588b6c26f4d5dc2a40364655391214441
+size 176201220
diff --git a/RU/bt2rad.m b/RU/bt2rad.m
new file mode 100644
index 0000000000000000000000000000000000000000..6cba838fae5b3e93754ca4bf268b60c8a4766682
--- /dev/null
+++ b/RU/bt2rad.m
@@ -0,0 +1,53 @@
+function radiance = bt2rad(freq,bt);
+%
+% function radiance = bt2rad(freq,bt);
+%
+% compute radiance given wavenumbers and brightness temperature.
+%
+% Inputs:
+% freq wavenumbers (Nx1) in 1/cm
+% bt brightnes temperature (Nx1) in Kelvin
+% Outputs:
+% radiance Planck radiances (Nx1) in mW/(m^2 sr. 1/cm)
+%
+% see also RADTOT.M, RAD2BT.M, TTORAD.M
+%
+% DCT 11/11-99
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+%
+
+% fundamental constants:
+% (Cohen, E.R. and B.N. Taylor, The 1986 CODATA recommended values
+% of the fundamental physical constants, Journal of Research of
+% the National Bureau of Standard, 92(2), March-April 1987.)
+
+h = 6.6260755E-34; % Planck constant in Js
+c = 2.99792458E8; % photon speed in m/s
+k = 1.380658E-23; % Boltzmann constant in J/K
+
+c1 = 2*h*c*c*1e8;
+c2 = h*c/k*1e2;
+
+radiance = 1e3 * c1.*(freq.*freq.*freq)./(exp((c2.*freq)./bt)-1);
+
+ind = find(bt < 0);
+if length(ind) ~= 0
+ radiance = nan*real(radiance);
+ disp('WARNING: negative inputs to BT2RAD')
+end
+ind = find(~isreal(bt));
+if length(ind) ~= 0
+ radiance = nan*real(radiance);
+ disp('WARNING: imaginary inputs to BT2RAD')
+end
+ind = find(isnan(bt));
+if length(ind) ~= 0
+ radiance = nan*real(radiance);
+ disp('WARNING: NaN inputs to BT2RAD')
+end
+ind = find(isinf(bt));
+if length(ind) ~= 0
+ radiance = nan*real(radiance);
+ disp('WARNING: Inf inputs to BT2RAD')
+end
+return
diff --git a/RU/convert_anc_to_structs.m b/RU/convert_anc_to_structs.m
new file mode 100644
index 0000000000000000000000000000000000000000..b91bce391c7c51030c490570b79d8d582d913caf
--- /dev/null
+++ b/RU/convert_anc_to_structs.m
@@ -0,0 +1,83 @@
+function [ICT_Param,RU,nonlin,polcorr_param] = convert_anc_to_structs(nc_struct);
+%
+% function [ICT_Param,RU,nonlin,polcorr_param] = convert_anc_to_structs(nc_struct);
+%
+% This function converts the contents from the CrIS L1b RU ancillary file
+% to ICT_Param, RU, nonlin, and polcorr structures for use in the RU calculations
+%
+% Input
+% nc_struct structure of vars read from CrIS L1b RU ancillary file
+%
+% Output
+% ICT_Param ICT ancillary parameters [struct]
+% RU Radiometric Uncertainty ancillary parameters [struct]
+% nonlin nonlinearity corretion ancillary parameters [struct]
+% polcorr_param polarization correction ancillary parameters [struct]
+%
+% JKT 2021-Aug-26
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ICT_Param.Earth_Temperature = nc_struct.ICT_Param_Earth_Temperature;
+ICT_Param.lw_ICT_wnum = nc_struct.wnum_lw;
+ICT_Param.lw_ICT_emiss = nc_struct.ICT_Param_lw_ICT_emiss;
+ICT_Param.lw_ICT_Baffle_emiss = nc_struct.ICT_Param_lw_ICT_Baffle_emiss;
+ICT_Param.lw_Housing_emiss = nc_struct.ICT_Param_lw_Housing_emiss;
+ICT_Param.lw_ScanBaffle_emiss = nc_struct.ICT_Param_lw_ScanBaffle_emiss;
+ICT_Param.lw_ScanMirror_emiss = nc_struct.ICT_Param_lw_ScanMirror_emiss;
+ICT_Param.lw_Earth_emiss = nc_struct.ICT_Param_lw_Earth_emiss;
+ICT_Param.mw_ICT_wnum = nc_struct.wnum_mw;
+ICT_Param.mw_ICT_emiss = nc_struct.ICT_Param_mw_ICT_emiss;
+ICT_Param.mw_ICT_Baffle_emiss = nc_struct.ICT_Param_mw_ICT_Baffle_emiss;
+ICT_Param.mw_Housing_emiss = nc_struct.ICT_Param_mw_Housing_emiss;
+ICT_Param.mw_ScanBaffle_emiss = nc_struct.ICT_Param_mw_ScanBaffle_emiss;
+ICT_Param.mw_ScanMirror_emiss = nc_struct.ICT_Param_mw_ScanMirror_emiss;
+ICT_Param.mw_Earth_emiss = nc_struct.ICT_Param_mw_Earth_emiss;
+ICT_Param.sw_ICT_wnum = nc_struct.wnum_sw;
+ICT_Param.sw_ICT_emiss = nc_struct.ICT_Param_sw_ICT_emiss;
+ICT_Param.sw_ICT_Baffle_emiss = nc_struct.ICT_Param_sw_ICT_Baffle_emiss;
+ICT_Param.sw_Housing_emiss = nc_struct.ICT_Param_sw_Housing_emiss;
+ICT_Param.sw_ScanBaffle_emiss = nc_struct.ICT_Param_sw_ScanBaffle_emiss;
+ICT_Param.sw_ScanMirror_emiss = nc_struct.ICT_Param_sw_ScanMirror_emiss;
+ICT_Param.sw_Earth_emiss = nc_struct.ICT_Param_sw_Earth_emiss;
+
+ICT_Param.View_Factor.Space = nc_struct.ICT_Param_A_Space;
+ICT_Param.View_Factor.BeamSplitterCold = nc_struct.ICT_Param_A_BeamSplitterCold;
+ICT_Param.View_Factor.BeamSplitterWarm = nc_struct.ICT_Param_A_BeamSplitterWarm;
+ICT_Param.View_Factor.OMAandFrame = nc_struct.ICT_Param_A_OMAandFrame;
+ICT_Param.View_Factor.ICTBaffle = nc_struct.ICT_Param_A_ICTBaffle;
+ICT_Param.View_Factor.ScanBaffle = nc_struct.ICT_Param_A_ScanBaffle;
+
+RU.dT_ict = nc_struct.RU_dT_ict;
+RU.de_ict = nc_struct.RU_de_ict;
+RU.dT_ict_ref_meas = nc_struct.RU_dT_ict_ref_meas;
+RU.dT_ict_ref_model = nc_struct.RU_dT_ict_ref_model;
+RU.dprpt = nc_struct.RU_dprpt;
+RU.dalpha = nc_struct.RU_dalpha;
+RU.lw_da2 = nc_struct.RU_lw_da2;
+RU.mw_da2 = nc_struct.RU_mw_da2;
+RU.sw_da2 = nc_struct.RU_sw_da2;
+
+nonlin.lw_a2 = nc_struct.nonlin_lw_a2;
+nonlin.mw_a2 = nc_struct.nonlin_mw_a2;
+nonlin.sw_a2 = nc_struct.nonlin_sw_a2;
+
+polcorr_param.delta_SSM = nc_struct.polcorr_delta_SSM;
+
+polcorr_param.lw_wn = repmat(nc_struct.wnum_lw(:),1,9);
+polcorr_param.lw_prpt = nc_struct.polcorr_lw_prpt;
+polcorr_param.lw_alph = nc_struct.polcorr_lw_alph;
+polcorr_param.lw_delta_cbb = nc_struct.polcorr_delta_cbb;
+polcorr_param.lw_delta_hbb = nc_struct.polcorr_delta_hbb;
+polcorr_param.mw_wn = repmat(nc_struct.wnum_mw(:),1,9);
+polcorr_param.mw_prpt = nc_struct.polcorr_mw_prpt;
+polcorr_param.mw_alph = nc_struct.polcorr_mw_alph;
+polcorr_param.mw_delta_cbb = nc_struct.polcorr_delta_cbb;
+polcorr_param.mw_delta_hbb = nc_struct.polcorr_delta_hbb;
+polcorr_param.sw_wn = repmat(nc_struct.wnum_sw,1,9);
+polcorr_param.sw_prpt = nc_struct.polcorr_sw_prpt;
+polcorr_param.sw_alph = nc_struct.polcorr_sw_alph;
+polcorr_param.sw_delta_cbb = nc_struct.polcorr_delta_cbb;
+polcorr_param.sw_delta_hbb = nc_struct.polcorr_delta_hbb;
+
+return
diff --git a/RU/cris_calbudget_1gran_nc_wrapper.m b/RU/cris_calbudget_1gran_nc_wrapper.m
new file mode 100644
index 0000000000000000000000000000000000000000..363e2da7d2c57c5a2714c9973726b5a712b3e564
--- /dev/null
+++ b/RU/cris_calbudget_1gran_nc_wrapper.m
@@ -0,0 +1,60 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% This is a wrapper script for testing the generation of the CrIS NASA L1b for a single L1b granule
+% and utilizes the netcdf4 formatted parameter files
+% The user must set the following (within USER INPUT section of script)
+%
+% sat_name satellite name ('J1' or 'SNPP')
+% l1b_fname L1b granule name (including path if not in same directory as code)
+% saveTotalRU_FLAG set to 1 if you want to save the RSS RU to a mat file
+%
+% The primary contributors to the CrIS radiometric uncertainty are
+% ICT temperature uncertainty
+% ICT emissivity uncertainty
+% uncertainty in the the other temperatures in the ICT radiometric model (modeled and measured)
+% uncertainty in the quadratic coefficient NLC coefficient (a2)
+% polarization correction uncertainty
+% (sensor polarization angle alpha, and combined SSM and sensor polarization prpt)
+%
+% RU term 1: e_ict ICT rad model, emitted radiance
+% RU term 2: T_ict ICT rad model, emitted radiance
+% RU term 3: T_SSM_scan_mirror_baffle ICT rad model, reflected radiance (measured T)
+% RU term 4: T_ScanBaffleCorrection ICT rad model, reflected radiance (modeled T)
+% RU term 5: T_ICT_baffle ICT rad model, reflected radiance (measured T)
+% RU term 6a: T_OMA, T_Frame ICT rad model, reflected radiance (measured T)
+% RU term 6b: T_BS_warm ICT rad model, reflected radiance (measured T)
+% RU term 7: T_BS_cold ICT rad model, reflected radiance (measured T)
+% RU term 8: T_Space ICT rad model, reflected radiance (modeled T)
+% RU term 9: prpt polarization correction
+% RU term 10: alpha polarization correction
+% RU term 11: a2 nonlinearity correction
+%
+% Refer to the NASA CrIS L1b RU document and ATBD for details on the ICT radiometric model
+% Each reflected term in the ICT radiometric model also has a view factor and emissivity, but the uncertainty
+% in those quantities produce negligible contributions to the CrIS radiometric uncertainty and are not perturbed for the RU estimate
+% J1 has an improved ICT emissivity, and the view factors of some of the reflected terms are set to zero as a result
+%
+% JKT 2020-Sep-13
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% START: USER INPUT
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+sat_name = 'SNPP'; % set to 'J1' or 'SNPP'
+l1b_fname = 'SNDR.SNPP.CRIS.20191031T0012.m06.g003.L1B.std.v03_08.G.201219125854.nc'; % set to file name of v3 L1b granule
+%sat_name = 'J1'; % set to 'J1' or 'SNPP'
+%l1b_fname = 'SNDR.J1.CRIS.20191031T0012.m06.g003.L1B.std.v03_08.G.201204050506.nc'; % set to file name of v3 L1b granule
+l1b_pname = '/Users/joet/git/cris_l1b_user/RU/TEST_INPUT'; % set to path to v3 L1b granule
+ru_code_pname = pwd; % set to path of the RU code
+
+saveTotalRU_FLAG = 1; % set to 1 if you want to save the RSS RU to a mat file
+pname_out = fullfile(pwd,'RU_OUT_NC'); % set to path for RU output files
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% END: USER INPUT
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+mkdir(pname_out)
+addpath(ru_code_pname);
+
+RU = cris_gran_RU_ncparam(sat_name,fullfile(l1b_pname,l1b_fname),saveTotalRU_FLAG,pname_out);
diff --git a/RU/cris_gran_RU_ncparam.m b/RU/cris_gran_RU_ncparam.m
new file mode 100644
index 0000000000000000000000000000000000000000..a6e9c490f6390f50a9a4c8bf553c57198bfe8382
--- /dev/null
+++ b/RU/cris_gran_RU_ncparam.m
@@ -0,0 +1,284 @@
+function RU = cris_gran_RU_ncparam(sat_name,l1b_fname,saveTotalRU_FLAG,pname_out);
+%
+% function RU = cris_gran_RU_ncparam(sat_name,l1b_fname,saveTotalRU_FLAG,pname_out);
+%
+% This is a function for generation of the CrIS NASA L1b for a single L1b granule
+% and utilizes the netcdf4 formatted parameter files
+%
+% Inputs
+% sat_name satellite name ('J1' or 'SNPP')
+% l1b_fname L1b granule name (including path if not in same directory as code)
+% saveTotalRU_FLAG set to 1 if you want to save the RSS RU to a mat file
+% pname_out output directory
+%
+% JKT 2021-Aug-26
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+% code will loop over bands
+bands{1} = 'lw';
+bands{2} = 'mw';
+bands{3} = 'sw';
+
+%%============================================================
+%% RU STEP 1a: Load ancillary data from static file (satellite dependent)
+%%============================================================
+switch sat_name
+ case 'J1'
+ % J1
+ jnk = rd_nc4('J1_L1B_ru_tool_param_file.nc');
+ case 'SNPP'
+ % SNPP
+ jnk = rd_nc4('SNPP_L1B_ru_tool_param_file.nc');
+end
+%%============================================================
+%% RU STEP 1b: separate .nc param file contents into ICT_Param, RU,
+%% nonlin, polcorr_param structures used in this code package
+%%============================================================
+[ICT_Param,RU,nonlin,polcorr_param] = convert_anc_to_structs(jnk);
+
+%%============================================================
+%% RU STEP 1c: copy data from anc file to expected var names for this code package
+%%============================================================
+dT_ict = RU.dT_ict;
+de_ict = RU.de_ict;
+dT_ict_ref_meas = RU.dT_ict_ref_meas;
+dT_ict_ref_model = RU.dT_ict_ref_model;
+dprpt = RU.dprpt;
+dalpha = RU.dalpha;
+da2.lw = RU.lw_da2;
+da2.mw = RU.mw_da2;
+da2.sw = RU.sw_da2;
+
+sensor_all.lw.a2_now = nonlin.lw_a2;
+sensor_all.mw.a2_now = nonlin.mw_a2;
+sensor_all.sw.a2_now = nonlin.sw_a2;
+
+%%%%%%%%%%
+%% create an opts structure
+%% not needed to calc RU estimate, but used here to record configuration options
+%%%%%%%%%%
+opts.extraFLAG = 1;
+opts.l1b_fname = l1b_fname;
+opts.bands = bands;
+opts.sat_name = sat_name;
+
+%%%%%%%%%%
+% create output file name
+% not required to calc RU estimate, but used here for optional output
+%%%%%%%%%%
+[~,fname_out,~] = fileparts(l1b_fname);
+fname_out = sprintf('%s_RU',fname_out);
+
+%%============================================================
+%% RU STEP 2: Load v3 L1b data
+%%============================================================
+l1b_dat = read_l1b_for_RU(l1b_fname);
+
+%%%%%%%%%%
+%% save opts params
+%% not needed to calc RU estimate, but used here to record configuration options
+%%%%%%%%%%
+opts.neonl = l1b_dat.aux.neon_wlen;
+save(fullfile(pname_out,sprintf('%sopts.mat',fname_out)), 'opts');
+
+%%%%%%%%%%
+% store lat and lon data to RU structure
+% not required to calc RU estimate, but are useful for plotting and statistical analysis
+%%%%%%%%%%
+%RU.lat = l1b_dat.lat;
+%RU.lon = l1b_dat.lon;
+
+%%%%%%%%%%
+%% sparsely populated sensor and user structs are used by this code package
+%%%%%%%%%%
+user_all.lw.band = 'LW';
+user_all.lw.v = l1b_dat.wnum_lw;
+user_all.lw.output_range = [user_all.lw.v(1) user_all.lw.v(end)];
+
+user_all.mw.band = 'MW';
+user_all.mw.v = l1b_dat.wnum_mw;
+user_all.mw.output_range = [user_all.mw.v(1) user_all.mw.v(end)];
+
+user_all.sw.band = 'SW';
+user_all.sw.v = l1b_dat.wnum_sw;
+user_all.sw.output_range = [user_all.sw.v(1) user_all.sw.v(end)];
+
+sensor_all.lw.band = 'LW';
+sensor_all.mw.band = 'MW';
+sensor_all.sw.band = 'SW';
+
+% generate RU, loop over detector bands
+%%%%%%%%%%%%%%%%%%%%
+for jj = 1:length(bands)
+
+ clear rad1 rad2 user sensor scTime sci eng geo ru_rad
+ user = user_all.(lower(bands{jj}));
+ sensor = sensor_all.(lower(bands{jj}));
+
+ v = l1b_dat.(sprintf('wnum_%s',bands{jj}));
+
+ %%============================================================
+ %% RU STEP 3: Remove polarization correction from L1b v3 radiance
+ %%============================================================
+ % rad_nopolcorr = radiance without nominal polarization correction
+ % rad_polcorr = radiance with nominal polarization correction
+ rad_polcorr = l1b_dat.(sprintf('rad_%s',bands{jj}));
+ [rad_nopolcorr,pol_params_nopolcorr] = pol_correctionRU_L1bonly(v, l1b_dat, sensor, user, ...
+ ICT_Param, polcorr_param, sat_name, 'prpt', -2); % prpt = prpt * (1+ru_value)
+
+ %% for RU calcs, rad1 should be the radiance prior to polarization correction
+ rad1 = rad_nopolcorr;
+ l1b_dat.(sprintf('rad_%s',bands{jj})) = rad1;
+ l1b_dat.(sprintf('rad_%s_pc',bands{jj})) = rad_polcorr;
+
+ %%============================================================
+ %% RU STEP 4: Assess ICT radiometric model uncertainties (RU terms 1, 2, 3, 4, 5, 6a/b, 7, 8)
+ %%============================================================
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % temperature and emissivity perturbations via cris_l1b_recal
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ % RU term 1: e_ict
+ %%%%%%%%%%%%%%%%%%%%
+ [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, lower(bands{jj}), ICT_Param, 'e_ICT', de_ict);
+ ru_rad.de_ict = rad2 - rad1;
+ ru_delta.de_ict = de_ict;
+
+ % RU term 2: T_ict
+ %%%%%%%%%%%%%%%%%%%%
+ [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, lower(bands{jj}), ICT_Param, 'T_ICT', dT_ict);
+ ru_rad.dT_ict = rad2 - rad1;
+ ru_delta.dT_ict = dT_ict;
+
+ % RU term 3: SSM_scan_mirror_baffle
+ %%%%%%%%%%%%%%%%%%%%
+ [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, lower(bands{jj}), ICT_Param, ...
+ 'SSM_scan_mirror_baffle', dT_ict_ref_meas);
+ ru_rad.dT_SSM_scan_mirror_baffle = rad2 - rad1;
+ ru_delta.dT_SSM_scan_mirror_baffle = dT_ict_ref_meas;
+
+ % RU term 4: T_ScanBaffleCorrection
+ %%%%%%%%%%%%%%%%%%%%
+ [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, lower(bands{jj}), ICT_Param, ...
+ 'T_ScanBaffleCorrection', dT_ict_ref_model);
+ ru_rad.dT_ScanBaffleCorrection = rad2 - rad1;
+ ru_delta.dT_ScanBaffleCorrection = dT_ict_ref_model;
+
+ % RU term 5: T_ICT_baffle
+ %%%%%%%%%%%%%%%%%%%%
+ [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, lower(bands{jj}), ICT_Param, ...
+ 'T_ICT_baffle', dT_ict_ref_meas);
+ ru_rad.dT_ICT_baffle = rad2 - rad1;
+ ru_delta.dT_ICT_baffle = dT_ict_ref_meas;
+
+ % RU term 6a: T_OMA, T_Frame
+ %%%%%%%%%%%%%%%%%%%%
+ [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, lower(bands{jj}), ICT_Param, ...
+ 'T_OMA', dT_ict_ref_meas);
+ ru_rad.dT_OMA = rad2 - rad1;
+ ru_delta.dT_OMA = dT_ict_ref_meas;
+
+ % RU term 6b: T_BS_warm %% We could potentially combine term 6a and 6b
+ %%%%%%%%%%%%%%%%%%%%
+ [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, lower(bands{jj}), ICT_Param, ...
+ 'T_BS_warm', dT_ict_ref_meas);
+ ru_rad.dT_BS_warm = rad2 - rad1;
+ ru_delta.dT_BS_warm = dT_ict_ref_meas;
+
+ % RU term 7: T_BS_cold
+ %%%%%%%%%%%%%%%%%%%%
+ [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, lower(bands{jj}), ICT_Param, ...
+ 'T_BS_cold', dT_ict_ref_meas);
+ ru_rad.dT_BS_cold = rad2 - rad1;
+ ru_delta.dT_BS_cold = dT_ict_ref_meas;
+
+ % RU term 8: T_Space (Earth / Space reflected out.eng.ICT_Param.Earth_Temperature)
+ %%%%%%%%%%%%%%%%%%%%
+ [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, lower(bands{jj}), ICT_Param, ...
+ 'T_Space', dT_ict_ref_model);
+ ru_rad.dT_Earth = rad2 - rad1;
+ ru_delta.dT_Earth = dT_ict_ref_model;
+
+ %%============================================================
+ %% RU STEP 5: Assess RU associated with polarization correction (RU terms 9, 10)
+ %%============================================================
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % RU associated with polarization correction
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ % RU term 9: polarization correction, prpt
+ %%%%%%%%%%%%%%%%%%%%
+ %% apply polarization correction with perturbed prpt
+ [rad3,pol_params_prpt_pert] = pol_correctionRU_L1bonly(v, l1b_dat, sensor, user, ICT_Param, ...
+ polcorr_param, sat_name, 'prpt', dprpt);
+ ru_rad.dprpt = rad_polcorr - rad3;
+ ru_delta.dprpt = dprpt;
+
+ % RU term 10: polarization correction, alpha
+ %%%%%%%%%%%%%%%%%%%%
+ %% apply polarization correction with perturbed alpha
+ [rad3,pol_params_prpt_pert] = pol_correctionRU_L1bonly(v, l1b_dat, sensor, user, ...
+ ICT_Param, polcorr_param, sat_name, 'alpha', dalpha);
+ ru_rad.dalpha = rad_polcorr - rad3;
+ ru_delta.dalpha = dalpha;
+
+ %%============================================================
+ %% RU STEP 6: Assess RU associated with nonlinearity correction (RU term 11)
+ %%============================================================
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % RU associated with nonlinearity correction
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ % RU term 11: nonlinearity correction, a2
+ %%%%%%%%%%%%%%%%%%%%
+ ru_rad.dnlc = estimate_nlc_ru_L1bonly(sensor,user,l1b_dat,ICT_Param,da2);
+ ru_delta.dnlc = da2.(bands{jj});
+
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ % RU subtotals and total
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ ru_rad.refl_meas = sqrt(ru_rad.dT_SSM_scan_mirror_baffle.^2 + ...
+ ru_rad.dT_ICT_baffle.^2 + ...
+ ru_rad.dT_OMA.^2 + ...
+ ru_rad.dT_BS_warm.^2 + ...
+ ru_rad.dT_BS_cold.^2);
+
+ ru_rad.refl_mod = sqrt(ru_rad.dT_ScanBaffleCorrection.^2 + ...
+ ru_rad.dT_Earth.^2);
+
+ ru_rad.pol = sqrt(ru_rad.dprpt.^2 + ...
+ ru_rad.dalpha.^2);
+
+ %%============================================================
+ %% RU STEP 7: Combine RU contributions via RSS
+ %%============================================================
+ ru_rad.total = sqrt(ru_rad.de_ict.^2 + ...
+ ru_rad.dT_ict.^2 + ...
+ ru_rad.dT_SSM_scan_mirror_baffle.^2 + ...
+ ru_rad.dT_ICT_baffle.^2 + ...
+ ru_rad.dT_OMA.^2 + ...
+ ru_rad.dT_BS_warm.^2 + ...
+ ru_rad.dT_BS_cold.^2 + ...
+ ru_rad.dT_ScanBaffleCorrection.^2 + ...
+ ru_rad.dT_Earth.^2 + ...
+ ru_rad.dprpt.^2 + ...
+ ru_rad.dalpha.^2 + ...
+ ru_rad.dnlc.^2);
+
+ RU.(bands{jj}).v = v;
+ RU.(bands{jj}).rad = rad1;
+ RU.(bands{jj}).rad_polcorr = rad_polcorr;
+ RU.(bands{jj}).rad_nopolcorr = rad_nopolcorr;
+ RU.(bands{jj}).ru_rad = ru_rad;
+ RU.(bands{jj}).ru_delta = ru_delta;
+
+ RU_out.(bands{jj}).v = v;
+ RU_out.(bands{jj}).ru_delta = ru_delta;
+ RU_out.(bands{jj}).ru_rad_rss = RU.(bands{jj}).ru_rad.total;
+ RU_out.(bands{jj}).rad = RU.(bands{jj}).rad_polcorr;
+end %%%%%%%%%%%%%%%%%%%% end loop over bands, jj %%%%%%%%%%%%%%%%%%%%
+if saveTotalRU_FLAG == 1
+ save(fullfile(pname_out,sprintf('%s.mat',fname_out)), 'RU_out');
+end
+return
diff --git a/RU/cris_l1b_recal_l1b_only.m b/RU/cris_l1b_recal_l1b_only.m
new file mode 100644
index 0000000000000000000000000000000000000000..a879a02403597271f54c86398f476e7ba9458e2a
--- /dev/null
+++ b/RU/cris_l1b_recal_l1b_only.m
@@ -0,0 +1,75 @@
+function [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, band, ICT_Param, ru_term, ru_value);
+%
+% function [rad2,B1,B2] = cris_l1b_recal_l1b_only(user, sensor, l1b_dat, band, ICT_Param, ru_term, ru_value);
+%
+% recalibrate radiance data using perturbed ICT radiometric model temperatures and/or
+% emissivities
+%
+% based on recal.m (DCT/ROK 12/10/99) and adapted for use with CrIS ICT radiometric model
+%
+% Inputs:
+% user sparse user structure
+% .band
+% .v
+% .output_range
+% sensor sparse sensor structure
+% .a2_now
+% .band
+% l1b_dat l1b data structure needed for RU estimate
+% ICT_param ICT radiometric model parameters
+% ru_term string indicating ru_term to perturb
+% ru_value size of perturbation
+%
+% Outputs:
+% rad2 recalibrated radiance spectra
+% B1 ICT rad model details for last obs in l1b_dat, no perturbation
+% B2 ICT rad model details for last obs in l1b_dat, with perturbation
+%
+% Theory:
+% rad1 = Re{(C_e - C_ds)./(C_ict - C_ds)}.*(L_ict1 - L_ds1) + L_ds1
+% rad2 = Re{(C_e - C_ds)./(C_ict - C_ds)}.*(L_ict2 - L_ds2) + L_ds2
+% = ((rad1-L_ds1)./(L_ict1-L_ds1)) .* (L_ict2-L_ds2) + L_ds2
+%
+% where:
+% C_e = complex measured spectra (Earth)
+% C_ds = complex measured spectra (Deep Space)
+% C_ict = complex measured spectra (Internal Cal Target)
+% L_ict1 = predicted ICT radiance (no T/emiss perturbations)
+% L_ds1 = predicted DS radiance (no T/emiss perturbations)
+% L_ict2 = predicted ICT radiance (with perturbed T/emiss)
+% L_ds2 = predicted DS radiance (with perturbed T/emiss)
+%
+% and
+% [L_ict1,~] = predictICTrad_L1b(user2, obs_time_tai93, sci1, eng1, geo);
+% [L_ict2,~] = predictICTrad_L1b(user2, obs_time_tai93, sci2, eng2, geo);
+% Use of "user2" structures is a trick to compute on the user grid
+% user2 = sensor;user2.v_onaxis = user.v;
+%
+%
+% assuming perfect space view: L_ds1 = L_ds2 = 0
+%
+% rad2 = (rad1 ./ (L_ict1) .* L_ict2
+%
+% 2020-08-19 JKT adapted for use with NASA CrIS L1b data files, no need for L1a data
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+
+fprintf(1,'ICT rad model perturbation: %s, %f\n',ru_term,ru_value);
+rad1 = l1b_dat.(sprintf('rad_%s',band));
+
+[nchan,nfov,nfor,nscan] = size(rad1);
+ind1 = find(user.v == user.output_range(1));
+ind2 = find(user.v == user.output_range(2));
+
+user2 = sensor;user2.v_onaxis = user.v;
+rad2 = nan(size(rad1));
+[ictRad1,B1] = predictICTradRU_L1bonly(user2, l1b_dat, ICT_Param, 'none', 0);
+[ictRad2,B2] = predictICTradRU_L1bonly(user2, l1b_dat, ICT_Param, ru_term, ru_value);
+
+ictRad1 = ictRad1(ind1:ind2,:,:);
+ictRad2 = ictRad2(ind1:ind2,:,:);
+
+for ii = 1:nfor
+ rad2(:,:,ii,:) = (squeeze(rad1(:,:,ii,:)) ./ ictRad1) .* ictRad2;
+end
+
+return
diff --git a/RU/estimate_nlc_ru_L1bonly.m b/RU/estimate_nlc_ru_L1bonly.m
new file mode 100644
index 0000000000000000000000000000000000000000..0430b5723da15fbdc940c4524cd4652ec5ebb88f
--- /dev/null
+++ b/RU/estimate_nlc_ru_L1bonly.m
@@ -0,0 +1,141 @@
+function nlc_ru = estimate_nlc_ru_L1bonly(sensor,user,l1b_dat,ICT_Param,da2);
+%
+% function [nlc_ru,rad_out,rad_out_new] = estimate_nlc_ru_L1bonly(sensor,user,l1b_dat,ICT_Param,da2);
+%
+% estimate radiometric uncertainty associated with nonlinearity correction (NLC) quadratic coefficient (a2) uncertainty
+%
+% Inputs
+% sensor sparse sensor structure
+% .a2_now
+% .band
+% user sparse user structure
+% .band
+% .v
+% .output_range
+% l1b_dat l1b data structure needed for RU estimate
+% ICT_param ICT radiometric model parameters
+% da2 a2 uncertainty value
+%
+%
+% Outputs
+% nlc_ru radiometric uncertainty associated with quadratic coefficient (a2) uncertainty
+%
+% 2020-08-19 JKT
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+
+PGA_GainALL.lw = l1b_dat.aux.pga_gain_lw;
+PGA_GainALL.mw = l1b_dat.aux.pga_gain_mw;
+PGA_GainALL.sw = l1b_dat.aux.pga_gain_sw;
+
+band = lower(user.band);
+ind1 = find(user.v == user.output_range(1));
+ind2 = find(user.v == user.output_range(2));
+v = user.v(ind1:ind2);
+rad_in = l1b_dat.(sprintf('rad_%s',band));
+[nchan,nfov,nfor,nscan] = size(rad_in);
+
+avgIT = double(l1b_dat.aux.(sprintf('ict_real_mean_%s',band)) + sqrt(-1)*l1b_dat.aux.(sprintf('ict_imag_mean_%s',band)));
+avgSP = double(l1b_dat.aux.(sprintf('space_real_mean_%s',band)) + sqrt(-1)*l1b_dat.aux.(sprintf('space_imag_mean_%s',band)));
+avgIT_Vdc = double(l1b_dat.aux.(sprintf('ict_vdc_%s',lower(user.band))));
+avgSP_Vdc = double(l1b_dat.aux.(sprintf('space_vdc_%s',lower(user.band))));
+ES_Vdc = double(l1b_dat.aux.(sprintf('earth_vdc_%s',lower(user.band))));
+PGA_gain = PGA_GainALL.(band);
+delta_a2 = da2.(band);
+user2 = sensor;user2.v_onaxis = user.v;
+indx = find(user2.v_onaxis >= user.output_range(1) & user2.v_onaxis <= user.output_range(2));
+[ICTrad,~] = predictICTradRU_L1bonly(user2, l1b_dat, ICT_Param, 'none', 0);
+ICTrad = ICTrad(indx,:,:);
+
+% make ICTrad same size as rad_in
+%%%%%%%%%%%%%%%%%%%%
+ICTrad2 = nan(nchan,nfov,1,nscan);
+ICTrad2(:,:,1,:) = ICTrad;
+ICTrad2 = repmat(ICTrad2,[1,1,nfor,1]);
+clear ICTrad
+
+% make avgIT, avgSP same size as rad_in
+% avgIT(:,:,1) and avgSP(:,:,1) are same interferometer sweep direction as even FOR ES
+% avgIT(:,:,2) and avgSP(:,:,2) are same interferometer sweep direction as odd FOR ES
+%%%%%%%%%%%%%%%%%%%%
+%% ICT %%
+jnk = zeros(nchan,nfov,nfor);
+%% interferometer sweep direction 0
+scan_sweep_dir_indx = find(l1b_dat.scan_sweep_dir == 0);
+cal_view_sweep_dir_indx = find(l1b_dat.aux.cal_view_sweep_dir == 0);
+jnk(:,:,scan_sweep_dir_indx) = repmat(avgIT(:,:,cal_view_sweep_dir_indx),[1,1,15]);
+%% interferometer sweep direction 1
+scan_sweep_dir_indx = find(l1b_dat.scan_sweep_dir == 1);
+cal_view_sweep_dir_indx = find(l1b_dat.aux.cal_view_sweep_dir == 1);
+jnk(:,:,scan_sweep_dir_indx) = repmat(avgIT(:,:,cal_view_sweep_dir_indx),[1,1,15]);
+avgIT2 = repmat(jnk,[1,1,1,nscan]);
+%% DS %%
+jnk = zeros(nchan,nfov,nfor);
+%% interferometer sweep direction 0
+scan_sweep_dir_indx = find(l1b_dat.scan_sweep_dir == 0);
+cal_view_sweep_dir_indx = find(l1b_dat.aux.cal_view_sweep_dir == 0);
+jnk(:,:,scan_sweep_dir_indx) = repmat(avgSP(:,:,cal_view_sweep_dir_indx),[1,1,15]);
+%% interferometer sweep direction 1
+scan_sweep_dir_indx = find(l1b_dat.scan_sweep_dir == 1);
+cal_view_sweep_dir_indx = find(l1b_dat.aux.cal_view_sweep_dir == 1);
+jnk(:,:,scan_sweep_dir_indx) = repmat(avgSP(:,:,cal_view_sweep_dir_indx),[1,1,15]);
+avgSP2 = repmat(jnk,[1,1,1,nscan]);
+clear jnk avgIT avgSP
+
+% make avgIT_Vdc, avgSP_Vdc the same size as ES_Vdc
+%%%%%%%%%%%%%%%%%%%%
+%% ICT %%
+jnk = zeros(nfov,nfor,1);
+%% interferometer sweep direction 0
+scan_sweep_dir_indx = find(l1b_dat.scan_sweep_dir == 0);
+cal_view_sweep_dir_indx = find(l1b_dat.aux.cal_view_sweep_dir == 0);
+jnk(:,scan_sweep_dir_indx,1) = repmat(avgIT_Vdc(:,cal_view_sweep_dir_indx),[1,15]);
+%% interferometer sweep direction 1
+scan_sweep_dir_indx = find(l1b_dat.scan_sweep_dir == 1);
+cal_view_sweep_dir_indx = find(l1b_dat.aux.cal_view_sweep_dir == 1);
+jnk(:,scan_sweep_dir_indx,1) = repmat(avgIT_Vdc(:,cal_view_sweep_dir_indx),[1,15]);
+avgIT_Vdc2 = repmat(jnk,[1,1,nscan]);
+%% DS %%
+jnk = zeros(nfov,nfor,1);
+%% interferometer sweep direction 0
+scan_sweep_dir_indx = find(l1b_dat.scan_sweep_dir == 0);
+cal_view_sweep_dir_indx = find(l1b_dat.aux.cal_view_sweep_dir == 0);
+jnk(:,scan_sweep_dir_indx,1) = repmat(avgSP_Vdc(:,cal_view_sweep_dir_indx),[1,15]);
+%% interferometer sweep direction 1
+scan_sweep_dir_indx = find(l1b_dat.scan_sweep_dir == 1);
+cal_view_sweep_dir_indx = find(l1b_dat.aux.cal_view_sweep_dir == 1);
+jnk(:,scan_sweep_dir_indx,1) = repmat(avgSP_Vdc(:,cal_view_sweep_dir_indx),[1,15]);
+avgSP_Vdc2 = repmat(jnk,[1,1,nscan]);
+clear jnk
+
+% estimate complex ES spectra via simple complex calibration equation, input radiance, and mean DS and ICT spectra
+%%%%%%%%%%%%%%%%%%%%
+cSpc = (rad_in ./ ICTrad2) .* (avgIT2 - avgSP2) + avgSP2;
+
+% estimate complex ES spectra with perturbed nonlinearity correction
+%%%%%%%%%%%%%%%%%%%%
+a2 = sensor.a2_now(:);
+a2_new = a2(:) + delta_a2(:);
+
+cSpc_new = zeros(size(cSpc));
+avgIT2_new = zeros(size(avgIT2));
+avgSP2_new = zeros(size(avgSP2));
+
+for iFOV = 1:nfov
+ for iES = 1:nfor
+ for iSCAN = 1:nscan
+ cSpc_new(:,iFOV,iES,iSCAN) = (1 + 2*a2_new(iFOV)*ES_Vdc(iFOV,iES,iSCAN))/...
+ (1 + 2*a2(iFOV)*ES_Vdc(iFOV,iES,iSCAN)) .* cSpc(:,iFOV,iES,iSCAN);
+ avgIT2_new(:,iFOV,iES,iSCAN) = (1 + 2*a2_new(iFOV)*avgIT_Vdc2(iFOV,iES,iSCAN))/...
+ (1 + 2*a2(iFOV)*avgIT_Vdc2(iFOV,iES,iSCAN)) .* avgIT2(:,iFOV,iES,iSCAN);
+ avgSP2_new(:,iFOV,iES,iSCAN) = (1 + 2*a2_new(iFOV)*avgSP_Vdc2(iFOV,iES,iSCAN))/...
+ (1 + 2*a2(iFOV)*avgSP_Vdc2(iFOV,iES,iSCAN)) .* avgSP2(:,iFOV,iES,iSCAN);
+ end
+ end
+end
+
+rad_out = ICTrad2 .* (cSpc - avgSP2) ./ (avgIT2 - avgSP2); % rad_out - rad_in should equal zero within noise
+rad_out_new = ICTrad2 .* (cSpc_new - avgSP2_new) ./ (avgIT2_new - avgSP2_new);
+
+nlc_ru = real(rad_out_new - rad_out);
+
+return;
diff --git a/RU/pol_correctionRU_L1bonly.m b/RU/pol_correctionRU_L1bonly.m
new file mode 100644
index 0000000000000000000000000000000000000000..a4aedcd882642d554ac10a5d40cd53bd34c956f0
--- /dev/null
+++ b/RU/pol_correctionRU_L1bonly.m
@@ -0,0 +1,121 @@
+function [L_polcorr,extra] = pol_correctionRU_L1bonly(v, l1b_dat, sensor, user, ICT_Param, polcorr_param, sat_name, ru_term, ru_value)
+%
+% function [L_polcorr,extra] = pol_correctionRU_L1bonly(v, l1b_dat, sensor, user, ICT_Param, polcorr_param, sat_name, ru_term, ru_value)
+%
+% Perform polarization correction, assuming the scene select mirror and sensor behave as
+% a linear polarizer-analyzer pair
+% apply perturbation to prpt or alpha if ru_term, ru_value are set
+%
+% Inputs
+% v wavenumber scale
+% l1b_dat uncorrected calibrated radiance
+% sensor reduced sensor parameters structure
+% user reduced user parameters structure
+% ICT_Param ICT parameters (view factor, emissivity, etc)
+% pol_corr_param ICT parameters (view factor, emissivity, etc)
+% sat_name 'J1' or 'SNPP'
+% ru_term term to perturb, 'prpt' or 'alpha' [optional]
+% ru_value size of perturbation [optional, not valid if ru_term is empty]
+%
+% Outputs
+% L_polcorr polarization correction calibrated radiance
+% extra polarization parameter info
+%
+% JKT 13-Sep-2021
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+%
+
+%% Initialization
+%% -------------------------------------------------------------------------
+
+if isempty(ru_term)
+ ru_term = 'none';
+ ru_value = 0;
+end
+
+L = l1b_dat.(sprintf('rad_%s',lower(sensor.band)));
+
+%SSMflag = 12;
+%prpt_approach_label = 'approach5';
+
+%% Set delta angles (FOV and FOR dependent)
+%% -------------------------------------------------------------------------
+delta_SSM = polcorr_param.delta_SSM; % already in radians
+
+%% calculate ICT radiance
+%% -------------------------------------------------------------------------
+user2 = sensor;
+user2.v_onaxis = user.v;
+ind1 = find(user.v == user.output_range(1));
+ind2 = find(user.v == user.output_range(2));
+user2.v_onaxis = user2.v_onaxis(ind1:ind2);
+[ictRad,~] = predictICTradRU_L1bonly(user2, l1b_dat, ICT_Param, 'none', 0);
+
+T_cbb = 4;
+L_cbb = bt2rad(user2.v_onaxis,T_cbb);
+
+L_polcorr = nan(size(L));
+alph = nan(9,1);
+prpt = nan(size(L,1),9);
+
+prpt_wn = polcorr_param.(sprintf('%s_wn',lower(sensor.band)));
+prpt_allFOV = polcorr_param.(sprintf('%s_prpt',lower(sensor.band)));
+alph = polcorr_param.(sprintf('%s_alph',lower(sensor.band)));
+delta_cbb = polcorr_param.(sprintf('%s_delta_cbb',lower(sensor.band)));
+delta_hbb = polcorr_param.(sprintf('%s_delta_hbb',lower(sensor.band)));
+
+%% Load polarization correction parameters (prpt, alpha, delta_cbb, delta_hbb)
+%% And spectrally smooth prpt if ksmooth ~= 0
+%% -------------------------------------------------------------------------
+for iFOV = 1:9
+
+ if ~isempty(find(isnan(prpt_allFOV)))
+ error('NaN values in prpt data');
+ end
+ % polarization prpt parameter is on NSR sampling, interpolate to current wavenumber grid if needed
+ prpt_in = interp1(prpt_wn(:,iFOV),prpt_allFOV(:,iFOV),v,'linear','extrap');
+ prpt(:,iFOV) = prpt_in;
+
+end
+
+%% Apply RU perturbation if needed
+%% -------------------------------------------------------------------------
+switch ru_term
+ case 'prpt'
+ prpt = prpt * (1+ru_value); % perturb prpt by a percentage defined by ru_value
+ fprintf(1,'prpt perturbation: prpt = prpt * %2.2f\n',1+ru_value);
+ case 'alpha'
+ alph = alph + ru_value; % perturb alpha by ru_value (radians)
+ fprintf(1,'alpha perturbation: +%1.3f rad\n',ru_value);
+ otherwise
+ fprintf(1,'nominal polarization correction\n');
+end
+
+%% Calculate polarization bias in calibrated radiance E_p
+%% -------------------------------------------------------------------------
+E_p = nan(size(L));
+for iSCAN = 1:size(L,4)
+ for iFOV = 1:9
+ L_ext = squeeze(L(:,iFOV,:,iSCAN));
+ %% J1 does not have SSM_scan_mirror temperature
+ %% so approximate SSM_scan_mirror temperature with OMA_structure_input temperature
+ if isempty(strfind(lower(sat_name),'j1'))
+ B_ssm = bt2rad(user2.v_onaxis,l1b_dat.aux.ssm_scan_mirror_temp(iSCAN));
+ else
+ B_ssm = bt2rad(user2.v_onaxis,l1b_dat.aux.oma_struct_temp_1(iSCAN));
+ end
+ E_p(:,iFOV,:,iSCAN) = Ep_approx(prpt(:,iFOV), alph(iFOV), delta_SSM(iFOV,:), L_cbb,...
+ squeeze(ictRad(:,iFOV,iSCAN)), B_ssm, L_ext, delta_cbb(iFOV), delta_hbb(iFOV));
+ end
+end
+
+%% Subtract polarization bias from uncorrected calibrated radiance
+%% -------------------------------------------------------------------------
+L_polcorr = L - E_p;
+
+extra.prpt = prpt;
+extra.v = v;
+extra.alpha = alph;
+extra.delta_SSM = delta_SSM;
+
+return
diff --git a/RU/predictICTradRU_L1bonly.m b/RU/predictICTradRU_L1bonly.m
new file mode 100644
index 0000000000000000000000000000000000000000..1852584e77cc9f8a98bbd32af30e924ce314bc6a
--- /dev/null
+++ b/RU/predictICTradRU_L1bonly.m
@@ -0,0 +1,57 @@
+function [ICTrad,B] = predictICTradRU_L1bonly(sensor,l1b_dat,ICT_Param,ru_term,ru_value);
+% function [ICTrad,B] = predictICTradRU_L1bonly(sensor,l1b_dat,ICT_Param,ru_term,ru_value);
+%
+%
+% Inputs
+% sensor sparse spectral param structure
+% .v
+% .band
+% .output_range
+% l1b_dat l1b data structure needed for RU estimate
+% ICT_param ICT radiometric model parameters
+% ru_term string indicating ru_term to perturb
+% ru_value size of perturbation
+%
+% Outputs
+% ICTrad nCHAN x nFOV array of predicted ICT radiance on the
+% B ICT rad model details for last obs in l1b_dat
+%
+% 2020-08-19 JKT adapted for use with NASA CrIS L1b data files, no need for L1a data
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+
+band = lower(sensor.band);
+
+% get key dimensions
+[~,nFOV,nFOR,nSCAN] = size(l1b_dat.rad_lw);
+nCHAN = length(sensor.v_onaxis);
+%nSCAN = length(l1b_dat.aux.ict_prt1_temp);
+%nFOR = size(l1b_dat.obs_time_tai93,1);
+
+% initialize output arrays
+ICTrad = ones(nCHAN, nFOV, nSCAN) * NaN;
+
+%% Loop over scans and compute ICT radiance
+for si = 1 : nSCAN % loop on scans
+
+ % check that this row has some ES's
+ if isnan(max(l1b_dat.obs_time_tai93(1:nFOR, si)));continue;end
+
+ T_ICT = (l1b_dat.aux.ict_prt1_temp(si) + l1b_dat.aux.ict_prt2_temp(si)) / 2;
+ T_CRIS.oma_struct_temp_1 = l1b_dat.aux.oma_struct_temp_1(si);
+ T_CRIS.oma_struct_temp_2 = l1b_dat.aux.oma_struct_temp_2(si);
+ T_CRIS.ssm_scan_mirror_temp = l1b_dat.aux.ssm_scan_mirror_temp(si);
+ T_CRIS.bs_temp = l1b_dat.aux.bs_temp(si);
+ T_CRIS.ssm_baffle_temp = l1b_dat.aux.ssm_baffle_temp(si);
+ T_CRIS.ssm_baffle_corr_temp = l1b_dat.aux.ssm_baffle_corr_temp(si);
+
+ % Compute ICT predicted radiance
+ B = ICTradModelRU_L1b_aux(sensor, T_ICT, T_CRIS, ICT_Param, ru_term, ru_value);
+
+ % Same values for all FOVs
+ for iFov = 1:nFOV;ICTrad(:,iFov,si) = B.total;end
+
+end
+B.ru_term = ru_term;
+B.ru_value = ru_value;
+
+return
diff --git a/RU/rd_nc4.m b/RU/rd_nc4.m
new file mode 100644
index 0000000000000000000000000000000000000000..d1140a3e10bb39375dd9659155095ef01152c882
--- /dev/null
+++ b/RU/rd_nc4.m
@@ -0,0 +1,43 @@
+function [root] = rd_nc4(fname);
+%
+% read data from netcdf4 file
+%
+% Inputs
+% fname file name with path if not in current directory
+%
+% Outputs
+% root structure containing contents of fname, organized by group
+%
+% JKT, University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+%
+
+root.info = ncinfo(fname);
+for ii = 1:length(root.info.Variables)
+ varname = root.info.Variables(ii).Name;
+ data = h5read(fname,['/' varname]);
+ root.(varname) = data;
+end
+
+for ii = 1:length(root.info.Groups)
+ group.name = root.info.Groups(ii).Name;
+ group.info = ncinfo(fname,group.name);
+ for vv = 1:length(group.info.Variables)
+ varname = group.info.Variables(vv).Name;
+ root.(group.name).(varname) = h5read(fname,['/' group.name '/' varname]);
+ end
+end
+
+
+for ii = 1:length(root.info.Groups)
+ if isfield(root.info.Groups(ii),'Groups')
+ for jj = 1:length(root.info.Groups(ii).Groups)
+ group.name = root.info.Groups(ii).Groups(jj).Name;
+ group.info = ncinfo(fname,['/' root.info.Groups(ii).Name '/' group.name]);
+ for vv = 1:length(group.info.Variables)
+ varname = group.info.Variables(vv).Name;
+ root.(root.info.Groups(ii).Name).(group.name).(varname) = h5read(fname,['/' root.info.Groups(ii).Name '/' group.name '/' varname]);
+ end
+ end
+ end
+end
+return
diff --git a/RU/read_l1b_for_RU.m b/RU/read_l1b_for_RU.m
new file mode 100644
index 0000000000000000000000000000000000000000..7d61b9a039cd5323c6b525535f488ca8e9bf9475
--- /dev/null
+++ b/RU/read_l1b_for_RU.m
@@ -0,0 +1,128 @@
+function d = read_l1b_for_RU(fname);
+%
+% function d = read_l1b_for_RU(fname);
+%
+% this function reads the variables needed for the RU estimate from a NASA CrIS L1b granule (v3)
+%
+% Inputs
+% fname filename with path
+%
+% Outputs
+% d data structure
+% .aux.measured_laser_wlen [scalar]
+% .aux.neon_wlen [scalar]
+% .aux.ict_prt1_temp [NSCAN x 1]
+% .aux.ict_prt2_temp [NSCAN x 1]
+% .aux.oma_struct_temp_1 [NSCAN x 1]
+% .aux.oma_struct_temp_2 [NSCAN x 1]
+% .aux.ssm_scan_mirror_temp [NSCAN x 1]
+% .aux.bs_temp [NSCAN x 1]
+% .aux.ssm_baffle_temp [NSCAN x 1]
+% .aux.ssm_baffle_corr_temp [NSCAN x 1]
+% .aux.pga_gain_lw [NFOV x 1]
+% .aux.pga_gain_mw [NFOV x 1]
+% .aux.pga_gain_sw [NFOV x 1]
+% .aux.ict_real_mean_lw [NCHAN_LW x NFOV x 2]
+% .aux.ict_real_mean_mw [NCHAN_MW x NFOV x 2]
+% .aux.ict_real_mean_sw [NCHAN_SW x NFOV x 2]
+% .aux.ict_imag_mean_lw [NCHAN_LW x NFOV x 2]
+% .aux.ict_imag_mean_mw [NCHAN_MW x NFOV x 2]
+% .aux.ict_imag_mean_sw [NCHAN_SW x NFOV x 2]
+% .aux.space_real_mean_lw [NCHAN_LW x NFOV x 2]
+% .aux.space_real_mean_mw [NCHAN_MW x NFOV x 2]
+% .aux.space_real_mean_sw [NCHAN_SW x NFOV x 2]
+% .aux.space_imag_mean_lw [NCHAN_LW x NFOV x 2]
+% .aux.space_imag_mean_mw [NCHAN_MW x NFOV x 2]
+% .aux.space_imag_mean_sw [NCHAN_SW x NFOV x 2]
+% .aux.earth_vdc_lw [NFOV x NFOR x NSCAN]
+% .aux.earth_vdc_mw [NFOV x NFOR x NSCAN]
+% .aux.earth_vdc_sw [NFOV x NFOR x NSCAN]
+% .aux.ict_vdc_lw [NFOV x 2]
+% .aux.ict_vdc_mw [NFOV x 2]
+% .aux.ict_vdc_sw [NFOV x 2]
+% .aux.space_vdc_lw [NFOV x 2]
+% .aux.space_vdc_mw [NFOV x 2]
+% .aux.space_vdc_sw [NFOV x 2]
+% .aux.cal_view_sweep_dir [2 x 1]
+% .wnum_lw [NCHAN_LW x 1]
+% .wnum_mw [NCHAN_MW x 1]
+% .wnum_sw [NCHAN_SW x 1]
+% .rad_lw [NCHAN_LW x NFOV x NFOR x NSCAN]
+% .rad_mw [NCHAN_MW x NFOV x NFOR x NSCAN]
+% .rad_sw [NCHAN_SW x NFOV x NFOR x NSCAN]
+% .obs_time_tai93 [NFOR x NSCAN]
+% .scan_sweep_dir [NFOR x 1]
+% .lat [NFOV x NFOR x NSCAN]
+% .lon [NFOV x NFOR x NSCAN]
+%
+% 2020-08-19 JKT
+% University of Wisconsin-Madison Space Science and Engineering Center (UW-SSEC)
+
+%% aux vars
+aux_var{1} = 'measured_laser_wlen'; % scalar
+aux_var{end+1} = 'neon_wlen'; % scalar
+aux_var{end+1} = 'ict_prt1_temp'; % NSCAN x 1
+aux_var{end+1} = 'ict_prt2_temp'; % NSCAN x 1
+aux_var{end+1} = 'oma_struct_temp_1'; % NSCAN x 1
+aux_var{end+1} = 'oma_struct_temp_2'; % NSCAN x 1
+aux_var{end+1} = 'ssm_scan_mirror_temp'; % NSCAN x 1
+aux_var{end+1} = 'bs_temp'; % NSCAN x 1
+aux_var{end+1} = 'ssm_baffle_temp'; % NSCAN x 1
+aux_var{end+1} = 'ssm_baffle_corr_temp'; % NSCAN x 1
+
+aux_var{end+1} = 'pga_gain_lw'; % NFOV x 1
+aux_var{end+1} = 'pga_gain_mw'; % NFOV x 1
+aux_var{end+1} = 'pga_gain_sw'; % NFOV x 1
+
+aux_var{end+1} = 'ict_real_mean_lw'; % NCHAN_LW x NFOV x 2
+aux_var{end+1} = 'ict_real_mean_mw'; % NCHAN_MW x NFOV x 2
+aux_var{end+1} = 'ict_real_mean_sw'; % NCHAN_SW x NFOV x 2
+aux_var{end+1} = 'ict_imag_mean_lw'; % NCHAN_LW x NFOV x 2
+aux_var{end+1} = 'ict_imag_mean_mw'; % NCHAN_MW x NFOV x 2
+aux_var{end+1} = 'ict_imag_mean_sw'; % NCHAN_SW x NFOV x 2
+
+aux_var{end+1} = 'space_real_mean_lw'; % NCHAN_LW x NFOV x 2
+aux_var{end+1} = 'space_real_mean_mw'; % NCHAN_MW x NFOV x 2
+aux_var{end+1} = 'space_real_mean_sw'; % NCHAN_SW x NFOV x 2
+aux_var{end+1} = 'space_imag_mean_lw'; % NCHAN_LW x NFOV x 2
+aux_var{end+1} = 'space_imag_mean_mw'; % NCHAN_MW x NFOV x 2
+aux_var{end+1} = 'space_imag_mean_sw'; % NCHAN_SW x NFOV x 2
+
+aux_var{end+1} = 'earth_vdc_lw'; % NFOV x NFOR x NSCAN
+aux_var{end+1} = 'earth_vdc_mw'; % NFOV x NFOR x NSCAN
+aux_var{end+1} = 'earth_vdc_sw'; % NFOV x NFOR x NSCAN
+
+aux_var{end+1} = 'ict_vdc_lw'; % NFOV x 2
+aux_var{end+1} = 'ict_vdc_mw'; % NFOV x 2
+aux_var{end+1} = 'ict_vdc_sw'; % NFOV x 2
+
+aux_var{end+1} = 'space_vdc_lw'; % NFOV x 2
+aux_var{end+1} = 'space_vdc_mw'; % NFOV x 2
+aux_var{end+1} = 'space_vdc_sw'; % NFOV x 2
+
+aux_var{end+1} = 'cal_view_sweep_dir'; % 2 x 1
+
+%% root group vars
+root_var{1} = 'wnum_lw'; % NCHAN_LW x 1
+root_var{end+1} = 'wnum_mw'; % NCHAN_MW x 1
+root_var{end+1} = 'wnum_sw'; % NCHAN_SW x 1
+
+root_var{end+1} = 'rad_lw'; % NCHAN_LW x NFOV x NFOR x NSCAN
+root_var{end+1} = 'rad_mw'; % NCHAN_MW x NFOV x NFOR x NSCAN
+root_var{end+1} = 'rad_sw'; % NCHAN_SW x NFOV x NFOR x NSCAN
+
+root_var{end+1} = 'obs_time_tai93'; % NFOR x NSCAN
+
+root_var{end+1} = 'scan_sweep_dir'; % NFOR x 1
+
+root_var{end+1} = 'lat'; % NFOV x NFOR x NSCAN
+root_var{end+1} = 'lon'; % NFOV x NFOR x NSCAN
+
+for ii = 1:length(aux_var)
+ d.aux.(aux_var{ii}) = ncread(fname,sprintf('aux/%s',aux_var{ii}));
+end
+for ii = 1:length(root_var)
+ d.(root_var{ii}) = ncread(fname,sprintf('%s',root_var{ii}));
+end
+
+return