#!/usr/bin/perl
my $usage = <<EOHDR;

    find_xtal.pl <DesiredFreq> [<Tolerance>]

    Find the crystal (and integral divisor(s)) to come within <Tolerance> Hz of 
    <DesiredFreq>.  Specify frequency as #.# or #m# for MHz, or #k# for kHz.
    If you don't specify <Tolerance>, 5% is assumed.

EOHDR
#
# 20130216 MCMason: Added tolerance, sort by error (asc) and funkiness within error.
#                   Find *all* divs within TOL range
# 20120429 MCMason: original version
#
#   TODO:   Trim report: If there are multiple ways to get same frequency, keep only
#           the top N of best one(s) (based on funkiness)
#   TODO:   Improve funkiness calculation: e.g. try to get something proportional
#           to # dividers & difficulty.  Perhaps something like $e * (log($f)/log(2))
#
use strict;
use warnings;
use autodie;

my $dbg_funkiness=0;

my $Freq = shift or die $usage;
$Freq = txt_2_MHz($Freq);
my $Tol = shift;
if (defined $Tol) {
    $Tol = txt_2_MHz($Tol);
}
else {
    $Tol = $Freq * .05;
}

my @crystals = read_crystals();
my @results;

sub compute {
    my ($Fx, $d) = @_;
    return undef unless $d;

    my $f = int($Fx / $d);
    my $err = abs($Freq - $f);
    return undef if $err > $Tol;

    my $ar = [ factorize($d) ];
    my $funk = factor_funkiness(@$ar);
    return [ $Fx, $d, $f, $err, $ar, $funk ];
}

for my $Fx (@crystals) {
    # 'perfect' divisor
    my $div = $Fx / $Freq;

    # Surrounding integral divisors
    my $d = int $div;
    while (my $ar = compute($Fx, $d)) {
        last if ! defined $ar;
        my ($xFx, $xd, $xf, $xerr, $xar, $xfunk) = @$ar;
        print "Fx:$xFx, d:$xd, f:$xf, err:$xerr, funk:$xfunk  (div:$div)\n" if $dbg_funkiness;
        push @results, $ar;
        --$d;
    }
    $d = int $div+1;
    while (my $ar = compute($Fx, $d)) {
        last if ! defined $ar;
        my ($xFx, $xd, $xf, $xerr, $xar, $xfunk) = @$ar;
        print "Fx:$xFx, d:$xd, f:$xf, err:$xerr, funk:$xfunk  (div:$div)\n" if $dbg_funkiness;
        push @results, $ar;
        ++$d;
    }
}

print "\n\n";
@results = sort { $$a[3] <=> $$b[3] or $$a[5] <=> $$b[5] } @results;

my $prev_freq = -1;
for my $ar (@results) {

    my $funk = $dbg_funkiness ? " \tfunkiness=$$ar[5]" : "";

    if ($prev_freq ne $$ar[2]) { 
        printf "%11.11s (%7s) = %11s / %4s  div = %-16s$funk\n",
                commify($$ar[2]), commify($$ar[3]), commify($$ar[0]), commify($$ar[1]),
                join(" * ", @{$$ar[4]});
        $prev_freq = $$ar[2];
    }
    else {
        printf "                      = %11s / %4s  div = %-16s$funk\n",
                commify($$ar[0]), commify($$ar[1]),
                join(" * ", factorize($$ar[1]));
    }
}

sub txt_2_MHz {
    my $t = shift;
    my $orig = $t;

    # Multiplier (default is MHz for /\d+\.\d*/)
    my $mult = 1000000;

    if ($t =~ /m/i) {
        #  1M8432 => 1,843,000
        $t =~ s/m/./;
        $mult = 1000000;
    }
    elsif ($t =~ /k/i) {
        # 32k768  =>    32,768
        $t =~ s/k/./;
        $mult = 1000;
    }

    if ($t =~ /\d+(\.\d*)?|\.\d+/) {
        return $t * $mult;
    }

    die "txt_2_MHz: Unexpected input '$orig' => '$t', can't determine frequency.";
}

sub read_crystals {
    open my $FH, '<', 'crystals.txt';
    my %freqs;
    while (<$FH>) {
        next if /^\s*#/;
        next if /^\s*$/;
        if (/([.0-9]+)/) {
            $freqs{$1*1_000_000}=0;
        }
    }
    return sort { $a <=> $b } keys %freqs;
}

sub commify {
    my $s = shift;
    $s =~ s/(\d)(\d\d\d)$/$1,$2/;
    $s =~ s/(\d)(\d\d\d,)/$1,$2/g;
    return $s;
}

sub factor_funkiness {
    my @factors = @_;
    my $funkiness = 0;
    print "funkiness(",join(", ", @factors),")\n" if $dbg_funkiness;
    for my $t (@factors) {
        my ($f, $e) = split /\^/, $t;
        $e=1 if ! defined $e;

        if ($f == 2) {
            $funkiness += .05*$e; # no change
        }
        elsif ($f == 3) {
            $funkiness += 0.5 * $e;
        }
        elsif ($f == 5) {
            $funkiness += 0.6 * $e;
        }
        elsif ($f < 32) {
            $funkiness += 1 * $e;
        }
        else {
            $funkiness += 2 * $e;
        }
        print "\t$t => $funkiness\n" if $dbg_funkiness;
    }

    print "\tfinal == $funkiness\n" if $dbg_funkiness;
    return $funkiness;
}

sub factorize {
    my $num = shift;
    my @factors = ();

    my $factor=2;
    my $exp=0;
    while ($num%$factor == 0) {
        ++$exp;
        $num /= $factor;
    }
    if ($exp > 1) { push @factors, "$factor^$exp" }
    elsif ($exp)  { push @factors, $factor }

    $factor=3;
    while ($factor*$factor <= $num) {
        $exp = 0;
        while ($num%$factor == 0) {
            ++$exp;
            $num /= $factor;
        }
        if ($exp > 1) { push @factors, "$factor^$exp" }
        elsif ($exp)  { push @factors, $factor }
        $factor += 2;
    }
    push @factors, $num if $num>1;

    return @factors;
}