# MATLAB code page

```% fixed constants
pday1=102970;
pday2=100955;
d=.7755; %must be in cm for correct conversion of E to V/m
rho=886;
g=9.8;
b=8.20E-3;

%% Drop2
eta2=1.858*10^-5;
V2=500.8;
trise2=[13.12 17.46 11.27 12.55 11.58];
tfall2=[18.17 19.46 19.86];
vrise2=(.5E-3)./trise2;
vfall2=(.5E-3)./tfall2;
vravg2=mean(vrise2)
vfavg2=mean(vfall2)

E2=(V2/(300*d))/(3.33E-5)
endterm2=(vravg2+vfavg2)/(E2*vfavg2);
a2=sqrt(((b/(2*pday1))^2)+ ((9*eta2*vfavg2)/(2*g*rho)))-(b/(2*pday1))
m2=(4/3)*pi*(a2^3)*rho

q2=g*m2*endterm2

% ERROR PROPAGATION
% standard deviation of mean for time data
sdmtr2=sqrt((sum((trise2-mean(trise2)).^2))/(length(trise2)*(length(trise2)-1)))
sdmtf2=sqrt((sum((tfall2-mean(tfall2)).^2))/(length(tfall2)*(length(tfall2)-1)))
% standard deviation of mean for velocities
sdmvr2=(vravg2/mean(trise2))*sdmtr2
sdmvf2=(vfavg2/mean(tfall2))*sdmtf2
sdma2=((9*eta2)/(2*g*rho))*((sdmvf2)/sqrt(((b/(2*pday1))^2)+((9*eta2*vfavg2)/(2*g*rho))))
% sdm of mass
sdmm2=4*pi*rho*(a2^2)*sdma2
% sdm of q
sdmq2=sqrt(((g*(vravg2+vfavg2)*sdmm2/(E2*vfavg2))^2)+((m2*g*sdmvr2/(E2*vfavg2))^2)...
+ ((m2*g*vravg2*sdmvf2/(E2*(vfavg2^2)))^2))

%% Drop9
eta9=1.855*10^-5;
V9=507;
trise9=[10.33, 9.05, 9.17];
tfall9=[10.14, 10.18, 10.05];
vrise9=(.5E-3)./trise9;
vfall9=(.5E-3)./tfall9;
vravg9=mean(vrise9)
vfavg9=mean(vfall9)

E9=(V9/(300*d))/(3.33E-5)
endterm9=(vravg9+vfavg9)/(E9*vfavg9);
a9=sqrt(((b/(2*pday2))^2)+ ((9*eta9*vfavg9)/(2*g*rho)))-(b/(2*pday2))
m9=(4/3)*pi*(a9^3)*rho

q9=g*m9*endterm9

% ERROR PROPAGATION
% standard deviation of mean for time data
sdmtr9=sqrt((sum((trise9-mean(trise9)).^2))/(length(trise9)*(length(trise9)-1)))
sdmtf9=sqrt((sum((tfall9-mean(tfall9)).^2))/(length(tfall9)*(length(tfall9)-1)))
% standard deviation of mean for velocities
sdmvr9=(vravg9/mean(trise9))*sdmtr9
sdmvf9=(vfavg9/mean(tfall9))*sdmtf9
sdma9=((9*eta9)/(2*g*rho))*((sdmvf9)/sqrt(((b/(2*pday2))^2)+ ((9*eta9*vfavg9)/(2*g*rho))))
% sdm of mass
sdmm9=4*pi*rho*(a9^2)*sdma9
% sdm of q
sdmq9=sqrt(((g*(vravg9+vfavg9)*sdmm9/(E9*vfavg9))^2)+((m9*g*sdmvr9/(E9*vfavg9))^2)...
+ ((m9*g*vravg9*sdmvf9/(E9*(vfavg9^2)))^2))

%% Drop11
eta11=1.861*10^-5;
V11=507;
atrise11=[22.67];
btrise11=[5.34, 5.70, 5.36];
tfall11=[5.31, 5.90, 5.93, 5.83];
avrise11=(.5E-3)./atrise11;
bvrise11=(.5E-3)./btrise11;
vfall11=(.5E-3)./tfall11;
avravg11=mean(avrise11)
bvravg11=mean(bvrise11)
vfavg11=mean(vfall11)

E11=(V11/(300*d))/(3.33E-5)
aendterm11=(avravg11+vfavg11)/(E11*vfavg11);
bendterm11=(bvravg11+vfavg11)/(E11*vfavg11);
a11=sqrt(((b/(2*pday2))^2)+ ((9*eta11*vfavg11)/(2*g*rho)))-(b/(2*pday2))
m11=(4/3)*pi*(a11^3)*rho

aq11=g*m11*aendterm11
bq11=g*m11*bendterm11

% ERROR PROPAGATION
% standard deviation of mean for time data
sdmatr11=0
sdmbtr11=sqrt((sum((btrise11-mean(btrise11)).^2))/(length(btrise11)*(length(btrise11)-1)))
sdmtf11=sqrt((sum((tfall11-mean(tfall11)).^2))/(length(tfall11)*(length(tfall11)-1)))
% standard deviation of mean for velocities
sdmavr11=0
sdmbvr11=(bvravg11/mean(btrise11))*sdmbtr11
sdmvf11=(vfavg11/mean(tfall11))*sdmtf11
sdma11=((9*eta11)/(2*g*rho))*((sdmvf11)/sqrt(((b/(2*pday2))^2)+ ((9*eta11*vfavg11)/(2*g*rho))))
% sdm of mass
sdmm11=4*pi*rho*(a11^2)*sdma11
% sdm of q
asdmq11=sqrt(((g*(avravg11+vfavg11)*sdmm11/(E11*vfavg11))^2)+((m11*g*sdmavr11/(E11*vfavg11))^2)...
+ ((m11*g*avravg11*sdmvf11/(E11*(vfavg11^2)))^2))
bsdmq11=sqrt(((g*(bvravg11+vfavg11)*sdmm11/(E11*vfavg11))^2)+((m11*g*sdmbvr11/(E11*vfavg11))^2)...
+ ((m11*g*bvravg11*sdmvf11/(E11*(vfavg11^2)))^2))

%% Drop12
eta12=1.862*10^-5;
V12=506.5;
atrise12=[8.93, 8.84];
btrise12=[3.48 3.68, 3.50, 3.65, 3.43];
ctrise12=[1.15, 1.15, 1.30, 1.21, 1.36];
tfall12=[22.00, 22.34, 25.62, 22.70, 20.45];
avrise12=(.5E-3)./atrise12;
bvrise12=(.5E-3)./btrise12;
cvrise12=(.5E-3)./ctrise12;
vfall12=(.5E-3)./tfall12;
avravg12=mean(avrise12)
bvravg12=mean(bvrise12)
cvravg12=mean(cvrise12)
vfavg12=mean(vfall12)

E12=(V12/(300*d))/(3.33E-5)
aendterm12=(avravg12+vfavg12)/(E12*vfavg12);
bendterm12=(bvravg12+vfavg12)/(E12*vfavg12);
cendterm12=(cvravg12+vfavg12)/(E12*vfavg12);
a12=sqrt(((b/(2*pday2))^2)+ ((9*eta12*vfavg12)/(2*g*rho)))-(b/(2*pday2))
m12=(4/3)*pi*(a12^3)*rho

aq12=g*m12*aendterm12
bq12=g*m12*bendterm12
cq12=g*m12*cendterm12

% ERROR PROPAGATION
% standard deviation of mean for time data
sdmatr12=sqrt((sum((atrise12-mean(atrise12)).^2))/(length(atrise12)*(length(atrise12)-1)))
sdmbtr12=sqrt((sum((btrise12-mean(btrise12)).^2))/(length(btrise12)*(length(btrise12)-1)))
sdmctr12=sqrt((sum((ctrise12-mean(ctrise12)).^2))/(length(ctrise12)*(length(ctrise12)-1)))
sdmtf12=sqrt((sum((tfall12-mean(tfall12)).^2))/(length(tfall12)*(length(tfall12)-1)))
% standard deviation of mean for velocities
sdmavr12=(avravg12/mean(atrise12))*sdmatr12
sdmbvr12=(bvravg12/mean(btrise12))*sdmbtr12
sdmcvr12=(cvravg12/mean(ctrise12))*sdmctr12
sdmvf12=(vfavg12/mean(tfall12))*sdmtf12
sdma12=((9*eta12)/(2*g*rho))*((sdmvf12)/sqrt(((b/(2*pday2))^2)+ ((9*eta12*vfavg12)/(2*g*rho))))
% sdm of mass
sdmm12=4*pi*rho*(a12^2)*sdma12
% sdm of q
asdmq12=sqrt(((g*(avravg12+vfavg12)*sdmm12/(E12*vfavg12))^2)+((m12*g*sdmavr12/(E12*vfavg12))^2)...
+ ((m12*g*avravg12*sdmvf12/(E12*(vfavg12^2)))^2))
bsdmq12=sqrt(((g*(bvravg12+vfavg12)*sdmm12/(E12*vfavg12))^2)+((m12*g*sdmbvr12/(E12*vfavg12))^2)...
+ ((m12*g*bvravg12*sdmvf12/(E12*(vfavg12^2)))^2))
csdmq12=sqrt(((g*(cvravg12+vfavg12)*sdmm12/(E12*vfavg12))^2)+((m12*g*sdmcvr12/(E12*vfavg12))^2)...
+ ((m12*g*cvravg12*sdmvf12/(E12*(vfavg12^2)))^2))

%% Drop14
eta14=1.862*10^-5;
V14=506.5;
atrise14=[1.81, 2.00, 1.58];
btrise14=[15.53, 20.33, 15.95, 17.81, 19.33];
tfall14=[7.08, 7.63, 7.75, 7.46, 7.87];
avrise14=(.5E-3)./atrise14;
bvrise14=(.5E-3)./btrise14;
vfall14=(.5E-3)./tfall14;
avravg14=mean(avrise14)
bvravg14=mean(bvrise14)
vfavg14=mean(vfall14)

E14=(V14/(300*d))/(3.33E-5)
aendterm14=(avravg14+vfavg14)/(E14*vfavg12);
bendterm14=(bvravg14+vfavg14)/(E14*vfavg12);
a14=sqrt(((b/(2*pday2))^2)+ ((9*eta14*vfavg14)/(2*g*rho)))-(b/(2*pday2))
m14=(4/3)*pi*(a14^3)*rho

aq14=g*m14*aendterm14
bq14=g*m14*bendterm14

% ERROR PROPAGATION
% standard deviation of mean for time data
sdmatr14=sqrt((sum((atrise14-mean(atrise14)).^2))/(length(atrise14)*(length(atrise14)-1)))
sdmbtr14=sqrt((sum((btrise14-mean(btrise14)).^2))/(length(btrise14)*(length(btrise14)-1)))
sdmtf14=sqrt((sum((tfall14-mean(tfall14)).^2))/(length(tfall14)*(length(tfall14)-1)))
% standard deviation of mean for velocities
sdmavr14=(avravg14/mean(atrise14))*sdmatr14
sdmbvr14=(bvravg14/mean(btrise14))*sdmbtr14
sdmvf14=(vfavg14/mean(tfall14))*sdmtf14
sdma14=((9*eta14)/(2*g*rho))*((sdmvf14)/sqrt(((b/(2*pday2))^2)+ ((9*eta14*vfavg14)/(2*g*rho))))
% sdm of mass
sdmm14=4*pi*rho*(a14^2)*sdma14
% sdm of q
asdmq14=sqrt(((g*(avravg14+vfavg14)*sdmm14/(E14*vfavg14))^2)+((m14*g*sdmavr14/(E14*vfavg14))^2)...
+ ((m14*g*avravg14*sdmvf14/(E14*(vfavg14^2)))^2))
bsdmq14=sqrt(((g*(bvravg14+vfavg14)*sdmm14/(E12*vfavg14))^2)+((m14*g*sdmbvr14/(E14*vfavg14))^2)...
+ ((m14*g*bvravg14*sdmvf14/(E14*(vfavg14^2)))^2))

%% Plots of time distributions

subplot(2,2,1), bar(trise2), title('Drop2 rise time')
subplot(2,2,2), bar(tfall11), title('Drop11 fall time')
subplot(2,2,3), bar(tfall12), title('Drop12 fall time')
subplot(2,2,4), bar(btrise14), title('Drop14 rise time')
```