prologues := 3;
outputformat := "svg";
defaultfont := "cmr10";
defaultscale := 1.1;
color Blue,Teal,Orange,Green,Grey,Light,Wash;
Blue := (0,84/255,159/255); Teal := (0,177/255,183/255);
Orange := (246/255,168/255,0); Green := (87/255,171/255,39/255);
Grey := (100/255,101/255,103/255); Light := (236/255,236/255,237/255);
Wash := (232/255,241/255,250/255);
% Format displayed numbers to hundredths using integer strings.
vardef numlabel(expr value)=
 save n,whole,frac; numeric n,whole,frac;
 n:=round(100abs value); whole:=floor(n/100); frac:=n mod 100;
 (if value<-.005: "-" & fi decimal whole
  if frac>0: & "." & if frac<10: "0" & fi decimal frac fi)
enddef;
pickup pencircle scaled .9;
outputtemplate := "MPT-14-006.svg";
input graph;
numeric n,t[],observed[],model[],low[],high[],valid[];
n:=0;
gdata("data/readings.dat", row,
 n:=n+1; t[n]:=scantokens row1; observed[n]:=scantokens row2;
 model[n]:=scantokens row3; low[n]:=scantokens row4;
 high[n]:=scantokens row5; valid[n]:=scantokens row6;
);
vardef mainpoint(expr x,y)=(60+40x,120+16y) enddef;
vardef residualpoint(expr x,y)=(60+40x,77+30y) enddef;
def interval_at(expr x,lo,hi)=
 draw mainpoint(x,lo)--mainpoint(x,hi) withcolor Grey;
 draw mainpoint(x,lo)+(-3,0)--mainpoint(x,lo)+(3,0) withcolor Grey;
 draw mainpoint(x,hi)+(-3,0)--mainpoint(x,hi)+(3,0) withcolor Grey;
enddef;
numeric count,total; count:=0; total:=0;
for j=1 upto n:
 if valid[j]=1: count:=count+1; total:=total+observed[j]-model[j]; fi
endfor;
def publication_axes=
 draw mainpoint(0,6)--mainpoint(0,0)--mainpoint(8,0) withcolor Grey;
 draw residualpoint(0,.8)--residualpoint(0,-.8)--residualpoint(8,-.8) withcolor Grey;
 draw residualpoint(0,0)--residualpoint(8,0) dashed evenly withcolor Grey;
 for x=0 upto 8:
  label.bot(decimal x,residualpoint(x,-.8));
 endfor;
 label.bot("Time (s)",(220,38));
 for y=0 step 2 until 6: label.lft(decimal y,mainpoint(0,y)); endfor;
 for y=-.8 step .8 until .8: label.lft(numlabel(y),residualpoint(0,y)); endfor;
 label.top("(a) Response (a.u.)",(220,213));
 label.top("(b) Residual (a.u.)",(220,100));
enddef;
beginfig(1);
fill unitsquare xscaled 420 yscaled 280 withcolor white;
% BEGIN DRAWING
publication_axes;
for j=1 upto n:
 if valid[j]=1:
  interval_at(t[j],low[j],high[j]);
  draw fullcircle scaled 5 shifted mainpoint(t[j],observed[j]) withcolor Blue;
  fill fullcircle scaled 4 shifted residualpoint(t[j],observed[j]-model[j])
   withcolor Teal;
 fi
endfor;
draw mainpoint(0,model[1])--mainpoint(8,model[n])
 dashed evenly withpen pencircle scaled 1.3 withcolor Orange;
label("Synthetic response and model" infont "cmbx10",(220,257));
draw fullcircle scaled 5 shifted (90,242) withcolor Blue;
label.rt("Observed + bounds",(98,242));
draw (263,242)--(283,242) dashed evenly withcolor Orange;
label.rt("Model",(291,242));
label("Synthetic data; missing at 3 s. Bounds are illustrative." infont "cmr8",(220,14));
message "MPT-14-006: valid=" & decimal count
 & "; residual sum=" & numlabel(total)
 & "; mean=" & decimal(total/count);
% END DRAWING
label.lft("MPT-14-006" infont "cmr8",(412,270));
setbounds currentpicture to unitsquare xscaled 420 yscaled 280;
endfig;
end.
