shebang.rb

require 'statsample'


def f2(n)
	return n.odd? ? (n * 3 + 1) / 2 : n / 2

end

def seq(n)
	l = [n]
	
	while (n != 1)
		n = f2(n)
		
		l << n
		
	end
	return l
end

def log2(x)
	return Math.log(x) / Math.log(2.0)
end

def count(ns)
	
	l = seq(n = ns.to_i(2))
	
	nw = ns.length

        cg = l.index { |x| x < l[0] }

        cm = (0..cg).max_by { |x| l[x] }
	c = l.size

	return {
		'ns' => ns,
		'ns2' => ns.reverse,
		'nw' => nw,
		'c' => c,
		'cg' => cg,
		'cm' => cm,
		'hc' => c.to_f / ns.length,
		'l' => l,
		'l0' => len0(ns),
		'l1' => len1(ns)
		}
end

def axes2(x2)
	return x2.inject({}) { |h, x| h.merge({x => 'axes x1y2' }) }
	
end

def x1y2()
	return $a2 if (!$a2.nil?)
	$a2 = axes2(['z', 'd', 'e', 'mp2', 'hc', 'pd', 'hc_d', 'hc_x', 'z2']) if ($a2.nil?)
	return $a2
end

def off()
	return {'ns' => nil, 'ns2' => nil}
end

def opts()
	return x1y2().merge(off())
end

def d1(s)
        c = s.split('').select { |x| x == '1' }.size
        return c
        
end

def midpt2(ns)
	w2 = d1(ns) / 2
	l = ns.split('')
	i = j = 0
	while (i < w2 && j < l.size)
		i += l[j].to_i
		j += 1
	end
	return j.to_f / ns.length
end

def midpt(ns)
	w2 = ns.length / 4
	l = ns.split('')
	i = j = 0
	while (i < w2 && j < l.size)
		i += l[j].to_i
		j += 1
	end
	return j.to_f / ns.length
end

def d(s)
        c = s.split('').select { |x| x == '1' }.size
	d = c.to_f / s.length
        return d
        
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 len1(ns)
	return len(ns, /0+/)
end

def len0(ns)
	return len(ns, /1+/)
end

def len01(ns)
	return len1(ns), len0(ns)
end

def e(ns)
	return len01(ns).flatten.size.to_f / ns.length
end


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 out(fno, a, a2 = {}, t = nil)

	$f = {} if ($f.nil?)
	if (!$f.member?(fno)) then
		$f[fno] = {'fn' => "out#{fno}.txt", 'c' => 0} 
		$f[fno]['f'] = File.open($f[fno]['fn'], 'w')
	end

	f = $f[fno]
	f['f'].puts(a.keys.join("\t")) if (f['c'] == 0)
	f['f'].puts(a.values.join("\t"))
	f['c'] += 1
	f['f'].flush
	
	return if ((f['c'] - 1) % 1000 != 0)
	
	d = (f['c'] / 2000).to_i
	d = 1 if (d == 0)
	
	f2 = File.open("gnuplot#{fno}.cmd", 'w')
	fn = f['fn']
	st = t.nil? ? 'unset title; ' : "set title '#{t}'; "
	k = "key top left box opaque"
	f2.puts("set ytics nomirror; set y2tics; set colors classic; set #{k}; #{st}; plot \\")
	c1, c2 = ":(column('w'))", 'linecolor palette' if (a.key?('w'))
	
	(a.keys - ['t']).each \
	{
		|x|
		next if (a2.member?(x) && a2[x].nil?)
		f2.puts("'#{fn}' using (column('t')):(column('#{x}'))#{c1} every #{d} with line lw 2 #{c2} #{a2[x]} title '#{x}', \\")
	}
	f2.puts
	f2.close
	
end

def fmt2(x)
	x = x.dup
	x.keys.each \
	{
		|k|
		x.delete(k) if (k.end_with?('s'))
	}
	
	return x
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 hide(a, x)
	return (a.member?(x) && a[x].nil?)
end

def outa(f, l, a = {}, t = '')
	
	k = outd(f, l)
	k1 = "key top left opaque"
	f.puts("set colors classic; set #{k1}; set title '#{t}'; ")
	f.puts("set ytics nomirror; set y2tics;")
	f.puts("plot \\")
	ct = ''
	k, ct = [k - ['t'], "(column('t')):"] if (k.member?('t') && !hide(a, 't'))
	k.each \
	{ 
		|x| 
		next if (hide(a, x))
		opt = a.fetch(x, '')
		opt += ' with line lw 2 ' if (!opt.include?('with') && !opt.include?('pt'))
		f.puts("'$dat' using #{ct}(column('#{x}')) #{opt} title '#{x}',\\") 
	}
	f.puts
	f.puts("set ytics")
	
#	f.puts("pause -1;")

end

def outafn(l, fn = nil, t = '', a = {})
	fn = 'gnuplot.cmd' if (fn.nil?)
	outa(f = File.open(fn, 'w'), l, a, t.sub("{", "\\{").sub("}", "\\}"))
	f.close
	$stderr.puts([fn, l.size, t].inspect)
end

def fmt1(l, f)

	w = l.max_by { |x| x.length }.length
	
	f.puts("set colors classic; set ytics nomirror; set y2tics; set key top right opaque;\\")
	f.puts("plot \\")
	f.puts("[0:#{l.size - 1}][0:#{w}] '-' matrix using 2:1:3 with image title '', \\")
	f.puts("'-' using (column('d')) with line axes x1y2 dt 2 lw 3, \\")
	f.puts("'-' using 1 with line title 'e' axes x1y2 dt 2 lw 3")
	l1 = []
	l.each_with_index \
	{
		|z, i|
		
		l1 << [d(z), e(z)]
		c = w - z.length
		z[z.length...w] = '0' * c
		
		f.puts(z.split('').join(' '))
	}
	f.puts("eof")
	f.puts("eof")
	f.puts("d")
	l1.each { |x| f.puts(x[0]) }
	f.puts("eof")
	l1.each { |x| f.puts(x[1]) }
	f.puts("eof")
	f.puts("set ytics;")
	
end

def grid(l, fn, t = '')
	
	f = File.open(fn, 'w')
	f.puts("set palette negative grayscale; unset colorbox; set title '#{t}';")
	
	fmt1(l.map { |x| x.reverse }, f)
	f.puts("set palette; set colorbox; set title;")
	f.close
	$stderr.puts([fn, l.size, t].inspect)

end

def fn(n2)
	return "gnuplot#{$n}-#{n2}.cmd"
end

def keys(l, ks)
	return l.map { |x| Hash[[ks, x.values_at(*ks)].transpose] }
end

def out1(h, k, a)
	l = h.values.flatten

	l = keys(l, l[0].keys - $ks - ['l', 'l0', 'l1'])

	l1 = l.map { |x| x['z'] }
	mn, mx = [l1.min, l1.max]
	
	l.each { |x| x['z2'] = (x['z'] - mn) / (mx - mn) }

	l.sort_by! { |x| x[k] }
	
	t = a.inspect
	
	outafn(l, fn(1), t, opts().merge({'t' => nil, 'z' => nil, 'z_d' => nil}))

	grid(l.map { |x| x['ns'] }, fn(2), t)
end

def insert(a, x)
	a['ks'].each \
	{ 
		|k| 
		
		a[k][x[k]] = a[k].fetch(x[k], []) + [x]
		
		mm = [a[k].keys.min, a[k].keys.max]
		
		a[k].delete_if { |k1,| !mm.member?(k1) }
	}
	
end

def delete(a, x)
	a['ks'].each \
	{
		|k|
		next if (!a[k].member?(x[k]))
		
		a[k][x[k]].delete(x)
		k1 = a[k].keys
		
		a[k].delete(x[k]) if (a[k][x[k]].empty?)
		k2 = a[k].keys
		
		raise [k1, k2].inspect if (k1.min != k2.min)
	}
end

def unique(a, x)
	a['ks'].each \
	{
		|k|
		return true if (a[k].fetch(x[k], []).size == 1) 
	}
	return false

end

def pass(x, h, j, bs, km, j1)
	
	if ($x.nil?) then
		$x = {}
		$x['ks'] = km.keys.select { |x| x.end_with?('_x') }.map { |x| x.chomp('_x') }
		$x['ks'].each { |x| $x[x] = {} }
	end
	
	j2 = (0...h[j].size).bsearch { |i| h[j][i]['z'] <= x['z']  }
	j2 = h[j].size if (j2.nil?)

	j1.replace(j2.to_s)
	h[j][j2, 0] = x
	insert($x, x)
	
	return if (h[j].size <= bs)
	
	i = h[j].size - 1
	i -= 1 while (unique($x, h[j][i]))
	delete($x, h[j][i])
	h[j][i, 1] = []
	
end

def add(bs, km, h, ns, l1, wt, t, wz, j1)

	return 0 if (ns.nil? || (a1 = count(ns)).nil?)

	j = ns.length
	h[j] = [] if (!h.member?(j))
	
	return 0 if (h[j].map { |x1| x1['ns'] }.member?(ns))

	x = adj(a1, wt, km).merge({'t' => t})

	if (!wz.empty?) then
		x.merge!(vars(x, wz))
		predict([x], 'cg', wz)
		x['cg_e'] = (x['cg'] - x['cg_y']).abs
	end

	pass(x, h, j, bs, km, j1)

	x1 = x.dup
	l1 << x1 if (!l1.nil?)
	
	out(1, wt.merge({'t' => x1['t']})) if (wt.member?('cg_e_a'))

	return 1
	
end

def initd(w)

	return(dense(w, rand))
end

def initdh(w)
	return(dense(w, 0.5))
end

def initb(w)
	r = rand
	return '1' + (1...w).map { rand < r ? '0' : '1' }.join
	
end

def initw(w)
	s = dense(w, 0.5)

	m = rand(w)
	i = rand(w - m) + 1
	
	s[i, m] = rand(2).to_s * m
	return s
end

def initr(w)
	
	s = '1'
	b = 0
	m = rand(w - 1) + 1
	while (s.length < w)
		
		z = w - s.length
		
		s += b.to_s * [rand(m - 1) + 1, z].min
		b ^= 1
		
	end
	return s
end

def init(w)
	m = rand(w)
	ns1 = dense(w, 0.5)
	ns2 = initr(w)

	ns1, ns2 = [ns2, ns1] if (rand(2) == 0)
	
	i = rand(w - m) + 1
	ns1[i, m] = ns2[i, m]
	
	return ns1
end

def initx(w)
	return '1' if (w == 1)
	ns = dense(w, 0.5)
	
	r1, r2 = (1..2).map { rand(w - 1) + 1 }.sort
	
	w2 = r2 - r1 + 1
	s1 = ['01', '10'][rand(2)]
	s = (s1 * ((w2 / 2) + 1))[0...w2]
	
	
	ns[r1..r2] = s
	return ns
end

def inity(w)

	return '1' * w if (w <= 2)
	w2 = w / 2
	w1 = rand(w2)
	
	ns = '1' + '0' * (w - w1 - 1) + '1' * w1
	
	l = (1..(w - w1 - 1)).to_a
	
	(w2 - w1 - 1).times { ns[l.delete_at(rand(l.size)), 1] = '1' }

	return ns
	
end

def mark(w, c)
	l = (0...w).to_a
	
	ns = '0' * w
	
	c.times { ns[l.delete_at(rand(l.size))] = '1' }
	return ns
end

def mark2(w, c2)
	w1 = rand(w2 = w / 2)
	
	ns1 = mark(w1, c = rand(w1 + 1))
	ns2 = mark(w - w1 - 1, c2 - c)
	
#	p([[ns1, ns1.length, c], [ns2, ns2.length, c2 - c]])
	return ns1 + ns2
end

def flip(ns)
	
	ns = ns.dup
	ns.length.times { |i| ns[i] = (ns[i].to_i ^ 1).to_s }
	return ns
	
end

def initm(w)
	return '1' * w if (w <= 2)
	ns = '1' + mark2(w, w / 2 - 1)
	return ns
end

def initmf(w)
	return '1' * w if (w <= 2)
	
	ns = '1' + flip(mark2(w, w / 2)).reverse
	return ns
end

def range(w)
	return [1] if (w == 1)
	a = (1..w).to_a

	l = (1..2).map { a.delete_at(rand(a.size)) }.sort
	
	return (l[0]..l[1]).to_a
end

def scatter(w)
	return [1] if (w == 1)
	a = (1..w).to_a

	c = rand(a.size) + 1
	l = (1..c).map { a.delete_at(rand(a.size)) }.sort
	
	return l
end

def mut(ns, l)
	ns2 = ns.dup
	
	l.each { |i| ns2[i, 1] = (ns[i, 1].to_i ^ 1).to_s }
	return ns2
end


def muta(ns)
	return mut(ns, range(ns.length - 1))
end

def mutb(ns)
	return mut(ns, scatter(ns.length - 1))
end

def remix(ns, l1)
	ns2 = ns.dup
	l2 = l1.shuffle
	
	l1.size.times { |i| ns2[l2[i], 1] = ns[l1[i], 1] }
	return ns2
end

def remixa(ns)
	return remix(ns, range(ns.length - 1))
end

def remixb(ns)
	return remix(ns, scatter(ns.length - 1))
end

def inner(ns1, ns2)
	a, b = [ns1, ns2].sort_by { |x| x.length }
	d = b.length - a.length
	i = rand(d + 1)
	
	c = b.dup
	c[i, a.length] = a
	
	ns1.replace(b)
	ns2.replace(c)
end

def mix(ns1, ns2)

	inner(ns1 = ns1.dup, ns2 = ns2.dup)
	return (0...ns1.length).map { |x| [ns1, ns2][rand(2)][x, 1] }.join
	
end

def cut(ns1, ns2)
	
	inner(ns1 = ns1.dup, ns2 = ns2.dup)
	i = rand(ns1.length)
	return ns1[0...i] + ns2[i..-1]
end

def ext(ns)

	ns = ns.dup
	$bd.times { ns[rand(ns.length) + 1, 0] = rand(2).to_s }
	return ns
end

def short(ns)
	return nil if (ns.length <= 50)
	ns = ns.dup
	$bd.times { ns[rand(ns.length - 1) + 1, 1] = '' }
	return ns
	
end

def stat(l)
	return 0, 1.0 if (l.empty?)
	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, sd
end

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

def avg(l)
	return l.empty? ? 0 : sum(l).to_f / l.size
end

def avg2(l)
	return l.empty? ? 0 : avg(l)
end

def weighted(ops, ar, ar2, ar1)

	ops.size.times \
	{
		|i|
		next if (ar1[i] == 0)
		ar[ops[i]] = (ar1[i].to_f / ar2[i]).round(3)
	}
	return rand(ops.size) if (ar.size != ops.size)

	ar.replace(Hash[ar.sort_by { |x| -x[1] }])

	wt = sum(ar.values)
	
	x = rand() * wt
	i = x1 = 0
	i, x1 = [i + 1, x1 + ar[ops[i]]] while (x1 < x)
	return i - 1
end

def top(l, c)
	
	return l[rand(l.size < c ? l.size : l.size / c)]
end

def linediffs(h, k)
	l = h.values.flatten
	
	kd = "#{k}_d"
	
	l.sort_by! { |x| x[k] }
	
	mn = l[0][k]
	mx = l[-1][k]
	
	l.size.times \
	{
		|i|
		y = mn + (i.to_f / (l.size - 1)) * (mx - mn)
		d = (y - l[i][k])
		l[i][kd] = d
	}
end

def dist(ks, h, wt)

	l = h.values.flatten
	ks.each \
	{
		|k|
		
		l1 = l.map { |x| x[k] }.compact
		next if (l1.empty?)
		
		a, sd = stat(l1)
		
		sd = 1.0 if (sd == 0)
		wt.merge!({"#{k}_a" => a, "#{k}_s" => sd})
	}
	
end

def absnorm(k, h, wt)
	ka = "#{k}_a"
	ks = "#{k}_s"
	kz = "#{k}_x"
	h.values.flatten.each \
	{
		|x|
		x[kz] = (x[k] - wt[ka]).abs / wt[ks]
	}
end

def adj(x, wt, km)
	z = 0.0
	
	km.keys.each \
	{ 
		|k| 
		if (!wt.member?("#{k}_a")) then
			z += x.fetch(k, 0)
			next
		end
		
		next if (!x.member?(k))
		
		z += (x[k] - wt["#{k}_a"]) / wt["#{k}_s"] * km[k] 
	}
	
	x['z'] = z
	raise [x['z'], wt].to_s if (!x['z'].finite?)
	return x
end

def rank(x)
	return -x.fetch('z', 0)
end

def bitdiff(ns1, ns2)
	mn = [ns1.length, ns2.length].min
	b1 = ns1[-mn..-1].to_i(2)
	b2 = ns2[-mn..-1].to_i(2)
	s = (b1 ^ b2).to_s(2)
	t = 0
	s.length.times { |i| t += s[i].to_i }
	return t
end

def prefixdiff(x, ns1, ns2, k, n = '')
	x["pd#{n}"] = [bitdiff(x[k], ns1), bitdiff(x[k], ns2)].min.to_f / [ns1.length + 1, ns2.length + 1].max
	
end

def prefixdiffs(h, k = 'ns', n = '')
	l = h.values.flatten.sort_by! { |x| x[k] }

	l.size.times \
	{
		|x|
		prefixdiff(l[x], x == 0 ? '' : l[x - 1][k], x == l.size - 1 ? '' : l[x + 1][k], k, n)
	}
end

def prefixdiffs2(h)
	prefixdiffs(h, 'ns', 1)
	prefixdiffs(h, 'ns2', 2)
	h.values.flatten.each { |x| x['pd'] = [x['pd1'], x['pd2']].min }
end


def auxvars(km, h, wt)
	
	prefixdiffs2(h) if (km.member?('pd'))

	l1 = []
	km2 = km.keys.select { |x| x.include?('_') }
	dist(km1 = km.keys - km2, h, wt)
	
	km2.each \
	{
		|k|
		if (k.end_with?('_d')) then
			k1 = k.chomp('_d')
			linediffs(h, k1)
		elsif (k.end_with?('_x')) then
			k1 = k.chomp('_x')
			dist([k1], h, wt) if (!km1.member?(k1))
			absnorm(k1, h, wt)
		end
	}
	
	dist(km2, h, wt)

	h.values.flatten.each { |x| adj(x, wt, km) }
	
	h.each { |k, l| l.sort_by! { |x| rank(x) } }
				
end

def try(h, ops, ar, ar1, ar2, w1, bc)
	
	begin
		r = weighted(ops, ar, ar2, ar1)
		
		opc = ops[r]

		op, p = opc.start_with?('init') ? [opc, 0] : [opc[0..-2], opc[-1..-1].to_i]

		ks = h.keys.sort_by { |x| rank(h[x][0]) }
		
	end while (p > 0 && ks.empty?)

	h.delete(ks.pop) if (ks.size > bc)
	
	case (p)
		when 0
			w = ks.empty? ? w1 : ks[rand(ks.size)]
			ns = send(op, w)
			
		when 1, 2
		
			arg = (1..p).map { top(ks, 5) }
			
			l = arg.map { |x| top(h[x], 5)['ns'] }
			ns = send(op, *l)
			ns[0, 1] = '1' if (!ns.nil?)
			
	end
	
#	raise op if (ns[0, 1] != '1')
	return ns, r
	
end

def log10(x)
	return x == 0 ? 0 : Math.log(x) / Math.log(10)
end

def lh(ns)
	nw2 = ns.length / 2
	return ns[nw2..-1], ns[0..nw2]
end

def dlh(ns)
	nw1, nw2 = lh(ns)
	d1 = d(nw1)
	d2 = d(nw2)
	return d(ns) - 0.5, d1 - 0.5, d2 - 0.5, d1 - d2
end

def elh(ns)
	nw1, nw2 = lh(ns)
	e1 = e(nw1)
	e2 = e(nw2)
	return e(ns) - 0.5, e1 - 0.5, e2 - 0.5, e1 - e2
end

def d4(ns)
	return d(ns[0..(ns.length / 4)])
end

def e4(ns)
	return e(ns[0..(ns.length / 4)])
end

def not0(x)
	return x == 0 ? 1 : x
end

def has(ks, p, ks1)

	p.replace({})
	
	ks1.each \
	{
		|k|
		p[k] = ks.nil? || ks.keys.index { |x| x.start_with?(k) }
	}
end

def data1(x, ks = nil)
	
	ns = x['ns']
	
	has(ks, p = {}, ['d', 's', 'ta', 'mp1', 'mp2', 'a', 'mx', 'd', 'e'])
	
	d, d1, d2, d12 = dlh(ns) if (p['d'] || p['s'] || p['ta'])
	
	mp = midpt(ns) - 0.5 if (p['s'] || p['mp1'])
	mp2 = midpt2(ns) - 0.5 if (p['ta'] || p['mp2'])

	a = {}
	
	if (p['s']) then
		s1 = mp.abs / not0(d12.abs)
		s3 = mp.abs * d12.abs
		s5 = -mp * d12

		a.merge!({
			's1' => s1,
			's2' => 1 / not0(s1),
			's4' => Math.log(1 / not0(s3)),
			's3' => s3,
			's5' => -s5,
			})
	end

	if (p['ta']) then
		
		t1 = mp2.abs / not0(d12.abs)
		t3 = mp2.abs * d12.abs
		t5 = -mp2 * d12
		
		a.merge!({
			'ta1' => t1,
			'ta2' => 1 / not0(t1),
			'ta4' => Math.log(1 / not0(t3)),
			'ta3' => t3,
			'ta5' => -t5,
			})
	end
	
	a1 = {}
	if (p['a']) then
		a0 = avg(x['l0'])
		a1 = avg(x['l1'])
		a01 = avg(x['l0'] + x['l1'])
		mx = log2(log2(ns.to_i(2).to_f))

		a1 = {
			'a0' => a0,
			'a1' => a1,
			'a01' => a01,
			'a01t' => (a0 + a1) / 2,
			'a01d' => a0 + a1 - 2 * a01,
			'a0m' => a0 / mx,
			'a1m' => a1 / mx,
			}
	end
	
	a2 = {}
	if (p['mx']) then
		mx0 = x['l0'].max
		mx0 = 0 if (mx0.nil?)
		mx1 = x['l1'].max
		mx01 = [mx0, mx1].max

		a2 = {
			'mx0' => mx0,
			'mx1' => mx1,
			'mx01' => mx01,
			}
	end

	a3 = {}
	if (p['d']) then
		a3 = {
			'd' => d,
			'da' => d.abs,
			'd1' => d1,
			'd1a' => d1.abs,
			'd2' => d2,
			'd2a' => d2.abs,
			'd12' => d12, 
			'd4' => d4(ns) - 0.5, 
			}
	end
	
	a4 = {}
	if (p['mp1']) then
		a4.merge!({
			'mp1' => mp,
			'mp1a' => mp.abs,
			})
	end
	if (p['mp2']) then
		a4.merge!({
			'mp2' => mp2,
			'mp2a' => mp2.abs,
			})
	end

	a5 = {}
	if (p['e']) then
		e, e1, e2, e12 = elh(ns)
		a5 = {
			'e' => e,
			'ea' => e.abs,
			'e1' => e1,
			'e1a' => e1.abs,
			'e2' => e2,
			'e2a' => e2.abs,
			'e12' => e12,
			'e4' => e4(ns) - 0.5,
			}
	end
	
	return a.merge(a1).merge(a2).merge(a3).merge(a4).merge(a5)
end

def stat2(k, l)
	a, sd = stat(l)
	return {"a#{k}" => a, "s#{k}" => sd} 
end

def runs(ns)

	ns = ns.reverse
	
	s = ''
	l = []
	ns.length.times \
	{
		|i|
		s += ns[i]
		l, s = [l + [[s, i]], ''] if (ns[i + 1] != ns[i])
	}
	
	return l
	
end

def filter(l, b = nil, t = nil, t2 = nil)
	l = l.select { |x| x[0][0] == b.to_s } if (!b.nil?)

	return l.sort_by { |x| [x[0].length, x[1]] }.reverse[0...t] if (!t.nil?)
	return l.select { |x| x[0].length >= t2 } if (!t2.nil?)
	return l
end

def diffs(l)
	a = stat2(nil, l)
	l.sort!
	m = l[l.size / 2]
	ma = avg2(l.select { |x| x > m })
	am = m.nil? ? 0 : a['a'] / not0(m)
	l2 = (1...l.size).map { |x| l[x] - l[x - 1] }
	a2 = stat2('2', l2)
	mx = l.max
	
	return {'m' => m.nil? ? 0 : m,
		'ma' => ma, 
		'am' => am,
		'mn' => l.empty? ? 0 : l.min, 
		'mx' => mx.nil? ? 0 : mx,
		'c' => l.size, 
	}.merge(a).merge(a2)

end

def number(a, w, b)
	Hash[[a.keys.map { |x| "#{x}-#{w}-#{b}" }, a.values].transpose]
end

def diffs1(ns, w, b = 1)

	l1 = runs(ns)

	case w
		when 1
			l1 = filter(l1, b, 3)
		when 2
			l1 = filter(l1, b, 10)
		when 3
			l1 = filter(l1, b, 15)
		when 4
			l1 = filter(l1, b, nil, 3)
	end
	
	return number(diffs(l1.map { |s, i| i.to_f / ns.length }), w, b)
end

def params()
	return (1..4).inject([]) { |t, w| t += [0, 1].map { |b| [w, b] }}
end

def vars(x, a)
	
	ks = a.keys.select { |x| x.include?('-') }
	ks1 = a.keys - ks
	
	a1 = data1(x, Hash[[ks1, a.values_at(*ks1)].transpose])
	
	ks.each \
	{
		|x|
		k, w, b = x.split('-')
		next if (!(w =~ /\d+/))
#		p([k, w, b])
		a1.merge!(diffs1(x, w.to_i, b))
	}
	return a1
		
end
	
def corr(l, y1, yp = nil)
	
	yp, ye = ["#{y1}_y", "#{y1}_e"] if (yp.nil?)
	xa = avg(l.map { |x| x[y1] })
	ya = avg(l.map { |x| x[yp] })

	tx = ty = txy = e = 0.0
	m = nil
	l.each \
	{
		|z|
		
		x = z[y1]
		y = z[yp]

		z[ye] = (x - y).abs if (!ye.nil?)
		txy += (x - xa) * (y - ya)
		tx += (x - xa) ** 2
		ty += (y - ya) ** 2
	}
	d = Math.sqrt(tx) * Math.sqrt(ty)
	$z = 1e-4
	r = txy / (d.abs < $z ? $z : d)
	return r
	
end

def dot(x, z)

	t = z['c']
	
	(z.keys - ['c']).each { |v| t += z[v] * x[v] }
	return t
		
end

def predict(l, vy, z)
	l.each { |x| x["#{vy}_y"] = dot(x, z) }

end

def fit(l, vx, vy)

	a = {}

	(vx + [vy]).each { |v| a[v] = l.map { |x| x[v] }.to_vector() }
	
	begin
		r = Statsample::Regression.multiple(a.to_dataset(), vy)
#	rescue Statsample::Regression::LinearDependency
#		return nil
	end
	
#	$stderr.puts(r.summary)

	z = r.coeffs.merge({'c' => r.constant})
	
	predict(l, vy, z)
	r = corr(l, vy)

	return r, z

end

def fit3(l, ks, k1, mn2 = 0.04, mn1 = 0.10)
	l2 = []
	ks.each \
	{
		|k|
		r = corr(l, k, k1)
		l2 << [k, r]
	}
	
	l1 = l2.select { |x| !x[1].finite? }
	raise l1.inspect if (!l1.empty?)

	l3 = []
	l2.sort_by { |x| -x[1].abs }.select { |x| x[1].abs > mn1 }.each \
	{
		|k, r|

		next if ([l3.map { |r1,| (r.abs - r1.abs).abs }.min, mn2].compact[0] < mn2)
		
		p([k, r])
		l3 << [r, k]
	}
	$stderr.puts("#{l3.size} / #{l2.size}")
	
	r, z = fit(l, l3.map { |x| x[1] }, k1)
	
	p([r, Hash[z.sort_by { |x| -x[1].abs }]])
	
	l1 = z.select { |k, v| !v.finite? }
	raise l1.inspect if (!l1.empty?)
	
	return z
end

def analyzer(h, wz)

	l = h.values.flatten

	k1 = l[0].keys
	
	l.each { |x| x.merge!(data1(x)) }
	params.each { |w, b| l.each { |x| x.merge!(diffs1(x['ns'], w, b)) }}
	$ks = l[0].keys - k1 if ($ks.nil?)
	wz.replace(fit3(l, $ks, 'cg'))
	
end

def opt(km, w1, l1, n, k)
	h = {}

	ops = ['init', 'initd', 'initb', 'initw', 'initr',  'initx', 
		'inity',
		'initm', 'initmf',
		'initdh',
		'muta1', 'mutb1', 'remixa1', 'remixb1',
		'mix2', 'cut2']
	ops += ['ext1', 'short1'] if (km.member?('w2'))

	bs = 1000
	bc = 20
	$bd = 10
	
	t = c1 = 0
	
	wt = {}
	
	ar1 = [0] * ops.size
	ar2 = ar1.dup
	ar = {}
	k0 = z0 = t1 = t2 = t3 = t4 = t21 = 0
	wz = {}
	
	
	loop \
	{
		ns, r = try(h, ops, ar, ar1, ar2, w1, bc)
		
		c = add(bs, km, h, ns, l1, wt, t, wz, j1 = '')
				
		
		c1 += c
		ar1[r] += c
		ar2[r] += 1
		
		if (c == 1 && c1 != 0 && (c1 % 250 == 0)) then
				
			auxvars(km, h, wt)
			
			analyzer(h, wz) if (c1 % 1e3.to_i == 0)
				
			t21 = log10(t) - log10(t1)
			out(2, {
				'hmn' => h.keys.min,
				'hmx' => h.keys.max,
				'hvs' => h.values.flatten.size,
				"#{k}0" => k0,
				'z0' => z0,
				'hc' => h.size,
				'j1' => j1.to_i,
				't1' => log10(t1),
				't2' => log10(t2),
				't3' => log10(t3),
				't4' => log10(t4),
				't21' => t21,
				't' => t}, 
				axes2(['z0', 't1', 't2', 't3', 't4', 't21']))
			out1(h, k, km.merge({'k' => k, 'c1' => c1})) if (c1 % 1e3.to_i == 0)
		end

		if (t != 0 && t % 1e3 == 0) then
			
			l = h.values.flatten
			x1 = l.max_by { |x| x.fetch('z', 0) }
			x2 = l.min_by { |x| x.fetch('z', 0) }
			
			t1 = x1['t']
			t2 = x2['t']
			t3 = (l2 = l.map { |x| x['t'] }).min
			t4 = avg(l2)
			
			z0 = x1['z']
			k0 = l.map { |x| x.fetch(k, 0) }.max
			
			msg = {'t' => t,
				't21' => (t - t1),
				'c1' => c1,
				'f' => (c1.to_f / t * 100).round(2),
				'hmn' => h.keys.min,
				'hmx' => h.keys.max,
				'hc' => h.size,
				}
			z = Hash[$x['ks'].map { |k1| [k1, Hash[$x[k1].map { |k2, v| [k2, v.size] }.sort]]}]
			$stderr.puts(msg)
			$stderr.puts("\t" + ar.inspect + "\n\t" + z.inspect) if (t % 2e4 == 0)
#			break if ((!n.nil? && t == n) || (t > 1e4 && t21 > 1.0))
		
		end
		
		t += 1
	}
	
	return h

end

def run(a, w = 200, c = 5e4.to_i, k = 'cg')
		
#	$n = 0 if ($n.nil?)
#	$n += 1
#	$f = nil
	h = opt(a, w, nil, nil, k)

end

a = {
	'hc_d' => 1,
	'hc_x' => 1,
	'cg_d' => 1,
	'cg_x' => 1,
	'cm_d' => 1,
	'cm_x' => 1,
	'pd' => 1,
	'cg_e' => 1
	}


run(a)