adapt6b.rb

def dense(w, d)
	w2 = w - 1
	
	a = (0...w2).to_a
	s = '0' * w2
	(1..(d * w - 1)).map { a.delete_at(rand(a.size)) }.each { |x| s[x, 1] = '1' }
	return ('1' + s)
	
end

def len(ns, p)
	l = ns.split(p)
	l = [] if (l.nil?)
	l.shift if (l[0] == '')
	return l.map { |x| x.length }
end

def len01(ns)
	return len(ns, /0+/), len(ns, /1+/)
end

def hist(l)
	a = {}
	l.each { |x| a[x] = a.fetch(x, 0) + 1 }
	return a
end

def hdiff(h1, h2)
	d = 0
	(h1.keys | h2.keys).each { |x| d += (h1.fetch(x, 0) - h2.fetch(x, 0)).abs }
	return d
end

def hdiff2(h2a, h2b)
	return hdiff(h2a[0], h2b[0]) + hdiff(h2a[1], h2b[1])
end

def outd(f, l)
	f.puts('$dat << eof')
	k = l[0].keys
	
	f.puts(k.join("\t"))
	l.each { |x| f.puts(x.values.join("\t")) }

	f.puts('eof')
	return k
end

def outa(f, l, a = {}, t = '')
	
	k = outd(f, l)
	f.puts("set colors classic; set key top right opaque; set title '#{t}'; ")
	f.puts("set ytics nomirror; set y2tics;")
	f.puts("plot \\")
	kt = k.member?('t')
	knw = k.member?('nw')
	kz = k.member?('z1')
	
	(k - ['t', 'nw']).each \
	{ 
		|x| 
		next if (a.member?(x) && a[x].nil?)
		
		cc = knw ? ":(column('nw'))" : nil
		ct = kt ? "(#{kz && x != 'z1' ? '-' : nil}column('t')):" : nil
		
		opt = ' with line lw 2 '
		opt += ' linecolor palette ' if (!cc.nil?)
		opt += a.fetch(x, '')
		f.puts("'$dat' using #{ct}(column('#{x}'))#{cc} title '#{x}' #{opt},\\") 
	}
	f.puts
#	f.puts("reset; pause -1;")

end

def outafn(l, a = {}, fn = nil, t = '')
	outa(f = File.open(fn.nil? ? 'gnuplot.cmd' : fn, 'w'), l, a, t)
	f.close
end

def hadd(a, b)
	(a.keys | b.keys).each { |x| a[x] = a.fetch(x, 0) + b.fetch(x, 0) }
end

def hdiv(a, d)
	a.keys.each { |x| a[x] = a[x].to_f / d }
end

def hist2(ns)
	l1, l0 = len01(ns)
	return hist(l0), hist(l1)
end

def hdiv2(h, c)
	hdiv(h = h.dup, c)
	return h
end

def hsmooth(c, nw)
	
	t0 = {}
	t1 = {}
	c.times \
	{
		a0, a1 = hist2(dense(nw, 0.5))

		hadd(t0, a0)
		hadd(t1, a1)
	}
	hdiv(t0, c)
	hdiv(t1, c)

	return t0, t1
end

def hsmooth2(c, nw, t01)
	t0, t1 = t01
	a0, a1 = hist2(dense(nw, 0.5))

	hadd(t0, a0)
	hadd(t1, a1)

	t0b = hdiv2(t0, c)
	t1b = hdiv2(t1, c)
	return t0b, t1b
end

def r()
	return ((Time.now.to_f - Time.now.to_i) * 1e6).to_i
end

def sum(l)
	l.inject { |t, x| t + x }.to_f
end

def stat(l)
	
	t = l.inject { |a, x| a + x }
	t2 = l.inject(0) { |a, x| a + (x ** 2) }
	
        c = l.size
        a = t.to_f / c
        z = t2.to_f / c - a ** 2
        sd = Math.sqrt(z < 0 ? 0 : z)
        raise if (sd.nan?)
	
        return {"a" => a, "s" => sd}
end

def runsd(a = nil, s1 = nil)

	return {'x' => 0, 'x2' => 0, 'l' => [], 'c' => 20} if (a.nil?) 
	
	a['l'] << s1
	a['x'] += s1
	a['x2'] += (s1 ** 2)
	
	if (a['l'].size == a['c'] + 1) then
		x1 = a['l'].shift
		a['x'] -= x1
		a['x2'] -= (x1 ** 2)
	end
	
	a1 = a['x'] / a['l'].size
	z = a['x2'] / a['l'].size - a1 ** 2
	sd = Math.sqrt(z)
	z1 = sd / s1
	return z1
	
end

def adapt(nw, t01, a1)
	
	l = []
	t = 1
	c = 20
	
	sd = runsd()
	
	t0 = t1 = nil
	while (t < 100) 
		t0, t1 = hsmooth2(t, nw, t01)
		l << {'t' => t, 's1' => (s1 = sum(t0.values + t1.values)), 'nw' => nw}


		z1 = runsd(sd, s1)
		
		l[-1]['z1'] = z1
		break if (t >= c && z1 < 0.005)
		t += 1
	end
	
	t01.replace([t0, t1])
	a1.merge!({'t' => nw, 's_1' => l[-1]['s1'], 't1' => t})
	return l
end

def adapt2(nw, t01, a1)
	
	l = []
	t = 1
	th = x = x2 = 0
	c = 20

	sd = runsd()
	
	while (t < 200)
		hd = hdiff2(hist2(dense(nw, 0.5)), t01)
		th += hd
		
		s1 = th / t
		l << {'t' => t, 's1' => s1, 'nw' => nw}
		
		z1 = runsd(sd, s1)
		
		l[-1]['z1'] = z1
		break if (t >= c && z1 < 0.005)
		t += 1
	end
	
	a1.merge!({'t' => nw, 's_2' => l[-1]['s1'], 't2' => t})
	return l
	
end

l = []
l2 = []
l1 = []

(1..40).each \
{
	|x|
	t01 = [{}, {}]
	l1 << {}
	
	l += adapt(nw = x * 20, t01, l1[-1])
	l << {}
	
	l2 += adapt2(nw, t01, l1[-1])
	l2 << {}
	
	$stderr.puts(l1[-1].inspect)
}

outafn(l, {'z1' => 'axes x1y2'})
outafn(l2, {'z1' => 'axes x1y2'}, 'gnuplot1.cmd')
outafn(l1, {'t1' => 'axes x1y2', 't2' => 'axes x1y2'}, 'gnuplot2.cmd')