function [wimedian,wiresidual,wimad] = mad_jro(jro_wi,npmad,n1,missing)
%%
%% This program is written for calculating the median, residual, and the
%% median absolute deviation from the JRO obs and the TIEGCM Wi at JRO
%% BA Emery, March 2011

wimedian = missing;
wiresidual = missing;
wimad = missing;

clear wi_npmad;

nmad = 0;
for i=n1:n1+npmad-1
 if (jro_wi(i) > missing) 
  nmad = nmad + 1;
  wi_npmad(nmad) = jro_wi(i);
 end
end
% Calculate median, residual, and mad
%if (nmad > 0)
% Want half the points before finding the median to avoid edges (esp in obs)
if (nmad > npmad/2)
 nmid = fix((nmad+1)/2);
%n1,nmad,wi_npmad
 wimedian = median(wi_npmad);
% y=mad(x,0 or 1) calculates the mean or median absolute deviations;  NaNs are treated as missing
% for a vector y, mad(x,1) = median(abs(x-median(x)))
 wiresidual = wi_npmad(nmid)-wimedian;
 wimad = mad(wi_npmad,1);
% Check assertion (passed!)
% wimad2 = median(abs(wi_npmad-wimedian));
% if (abs(wimad-wimad2) > 0.1)
%wimad,wimad2
% end
end


