+from site_cons.site_init import link_wtk_graph
+
+
+FIGURES = ['expt-sawtooth']
+
# Get the passed in environment.
Import('env')
-for gp in Glob('*.gp'):
- env.Gnuplot(gp)
+wtk_graph = link_wtk_graph(env)
+
+xml_fit = File('sawmodel.xml')
+expt_sawtooth_data = [
+ env.Command(
+ 'expt-Ls.dat',
+ [xml_fit],
+ r"sed -n 's/^ *<L>\(.*\)<\/L>/\1/p' $SOURCE > $TARGET"),
+ ]
+
+for fig in FIGURES:
+ asyfile = '%s.asy' % fig # static .asy file
+ pyfig = fig.replace('-', '_')
+ data = '%s_data' % (pyfig)
+ asydata = []
+ if data in globals(): # generated data dependencies
+ asydata = globals()[data]
+ env.Asymptote([asyfile, wtk_graph] + asydata)
# Pass back the modified environment.
Return('env')
--- /dev/null
+import wtk_graph;
+
+size(15cm,10cm,IgnoreAspect);
+
+scale(Linear, Linear);
+real k = infinity; /* spring constant in N/m, already accounted for (0.05) */
+real xscale=1e9;
+real fscale=1e12;
+
+real xmin = -5e-9;
+real xmax = 300e-9;
+real fmin = -100e-12;
+real fmax = 1000e-12;
+real fnmax = 300e-12;
+
+real kB = 1.3806504e-23;
+real T = 300;
+real p = 3.7e-10;
+real WLC(real x, real L, real p=p, real T=T, real kB=kB,
+ real xLmax=0.99, real fmax=fnmax) {
+ real f = infinity;
+ if (x < 0)
+ return infinity;
+ if (x/L <= xLmax)
+ f = kB*T/p*(0.25*(1/(1-x/L)^2 - 1) + x/L);
+ if (f > fnmax)
+ return infinity;
+ return f;
+}
+
+void graphSawtooth(string file="datafile", real k,
+ int xcol=0, int fcol=1,
+ real xscale=1, real fscale=1, real dx=0, real df=0,
+ pen p=red) {
+ file fin=input(file).line();
+ real[][] a=fin.dimension(0,0);
+ a=transpose(a);
+ real[] x=a[xcol];
+ real[] f=a[fcol];
+ x = x - f/k; /* Remove cantilever extension */
+ graphData(x=x, y=f, xscale=xscale, yscale=fscale,
+ dx=dx, dy=df, p=p, dots=false);
+}
+
+/* expt.dat already scaled to nm and pN */
+graphSawtooth("expt.dat", k=k, p=psoft);
+
+file fin = input('expt-Ls.dat').line();
+real[][] a = fin.dimension(0,0);
+a = transpose(a);
+real[] Ls = a[0];
+int i;
+for (i=0; i<Ls.length; ++i) {
+ real fn(real x) {
+ real f = WLC(x / xscale, L=Ls[i]);
+ if (f == infinity) {
+ write(x, fnmax*fscale * 1.1);
+ return fnmax*fscale * 1.1;
+ }
+ write(x, f*fscale);
+ return f*fscale;
+ }
+ draw(graph(fn, a=0, b=0.9*Ls[i]*xscale), p=phard);
+}
+
+xlimits(xmin*xscale, xmax*xscale, crop=true);
+ylimits(fmin*fscale, fmax*fscale, crop=true);
+
+xaxis(sLabel("Distance (nm)"), BottomTop, LeftTicks);
+yaxis(sLabel("Force (pN)"), LeftRight, RightTicks);
+label(sLabel("Poly-ubiquitin force curve"), point(N), N);
+++ /dev/null
-set terminal pdf font 'Times,7' size 9cm, 6cm
-set output 'fig.pdf'
-unset key
-set yrange [-100:1000]
-set xrange [-5:250]
-set samples 10000
-set xtics nomirror
-set ytics nomirror
-# border 3 = bottom & left
-set border 3
-set xlabel 'Extension (nm)'
-set ylabel 'Force (pN)'
-kB = 1.3806504e-23
-T = 300
-p = 3.7e-10
-f(x,L) = 1e12 * kB*T/p*(0.25*(1/(1-x/L)**2 - 1) + x/L)
-F(x,L) = (x>=0 && x < L && f(x,L) <= 300) ? f(x,L) : 1/0
-plot 'expt.d' with lines, \
- F(x,41.0063) lt 2, \
- F(x,65.2187) lt 2, \
- F(x,88.4682) lt 2, \
- F(x,110.442) lt 2, \
- F(x,133.067) lt 2, \
- F(x,155.498) lt 2, \
- F(x,177.802) lt 2, \
- F(x,199.098) lt 2, \
- F(x,221.89) lt 2
graphSawtooth("fig3.dat", k=k, xscale=xscale, fscale=fscale,
df=400e-12, p=phard);
-//xlimits(xmin, xmax);
-//ylimits(ymin, ymax);
xaxis(sLabel("Distance (nm)"), BottomTop, LeftTicks);
yaxis(sLabel("Force (pN)"), LeftRight, RightTicks);
-/*
-set xrange [-5:250]
-set yrange [-10:400]
-set ytics 0,100
-
-set multiplot
-
-set size 1,0.42;
-set origin 0.0,0.0;
-set xlabel 'Extension (nm)'
-plot 'fig1.d' using (($1-$2/0.05)*1e9):($2*1e12)
-
-set size 1,0.29;
-set origin 0.0,0.42;
-set bmargin 0
-unset xlabel
-set xtics format ""
-set ylabel 'Force (pN)'
-plot 'fig2.d' using (($1-$2/0.05)*1e9):($2*1e12)
-
-set origin 0.0,0.71;
-set bmargin 0
-set xtics format ""
-unset ylabel
-plot 'fig3.d' using (($1-$2/0.05)*1e9):($2*1e12)
-*/
+label(sLabel("Simulated force curves"), point(N), N);