#!/usr/bin/env perl

use warnings;
use strict;
use DateTime;

#
# Takes a sequences of DCF77 data bits (one minuite per line)
# and decrypts and decodes the meteotime data in them.
#
# Eg:
#
# 001001100011001001001100011011001000000010001101000001100000
# 001101101100110001001010011011001000000010001101000001100000
# 000000101011111001001110011001001000000010001101000001100000
#
# gives
#
# 2018-05-10 09:32 line        0 bits: 0100 0111 1010 0100 1111 1001
# region: FR Grenoble, cycle: Lows day 2, when: 11/05/2018
#   day        = slightly cloudy
#   night      = high altitude fog
#   temp       = 9⁰C
#   anomaly    = 0
#   wind speed = 3-4 Bft
#   wind dir   = SW
#
#
#
# no parity or continuity checks are done on the time data,
# we rely on the 16 bits of check data in the meteotime stream
# to detect errors.
#
#

######################################################### BITS

# print bits in groups of 6
sub b6($) {
    my $v = shift;
    return join( ' ', ( unpack "a6" x ( ( length($v) / 6 ) - 1 ) . "a*", $v ) );
}

# print bits in groups of 4
sub b4($) {
    my $v = shift;
    return join( ' ', ( unpack "a4" x ( ( length($v) / 4 ) - 1 ) . "a*", $v ) );
}

# big endian bits to integer
sub be_bin_to_int($) {
    return oct( '0b' . shift );
}

# little endian bits to integer
sub le_bin_to_int($) {
    return oct( '0b' . reverse(shift) );
}

# litte endian bcd bits to integer
sub le_bcd_to_int($) {
    my $v   = shift;
    my $ret = 0;
    my $c   = 1;
    for ( my $i = 0 ; $i < length($v) ; $i += 4 ) {
        $ret += $c * le_bin_to_int( substr( $v, $i, 4 ) );
        $c *= 10;
    }
    return $ret;
}

######################################################### UTILS

# extract ciphertext and key from three consecutive bit strings
sub extract($$$) {
    my ( $d0, $d1, $d2 ) = @_;
    
    my $data = "";
    my $key  = "";
    my $cest = "";

    $data .= substr( $d0, 2, 6 );
    $data .= substr( $d0, 9, 6 );
    $data .= substr( $d1, 1, 14 );
    $data .= substr( $d2, 1, 14 );

    $key .= substr( $d1, 21, 7 ) . '0';     # min
    $key .= substr( $d1, 29, 6 ) . '00';    # hour
    $key .= substr( $d1, 36, 6 ) . '00';    # day of month
    $key .= substr( $d1, 45, 5 );           # month
    $key .= substr( $d1, 42, 3 );           # day of week
    $key .= substr( $d1, 50, 8 );           # year

    $cest = substr( $d1, 17, 1 );           # cest

    return ( $data, $key, $cest );
}

######################################################### CRYPTO

# So this is almost DES with the following changes:
#
# 1) There's no initial or final permutation
# 2) the block size is 40 bits, half blocks are 20 bits
# 3) There's no permuted choice 1, permuted choice 2 is different.
# 4) Rounds 8,7 and 3 have a double rotation for the key
# 5) The expansion permutation is replaced with a choice followed by catenation.
# 6) There's an extra trvial permutation after the xor with subkey before the sboxes
# 7) There are only 5 sboxes, and they're different from DES.
# 8) The permutation is different.

# figuring all this out was rather miserable.

# ror a bitstring of 20 bits
sub ror($) {
    my $v = shift;

    return substr( $v, 19, 1 ) . substr( $v, 0, 19 );
}

# not quite DES lookup table for PC2
# 0-19 from Left Key
# 20-39 from Right Key
my @compress_lut = (
    31, 1,  18, 16, 28, 19, 9,  7,  13, 4,  39, 2,  21, 5, 0, 38,
    15, 37, 25, 3,  23, 33, 35, 29, 12, 22, 24, 34, 14, 10
);

# make 30 bit subkey from 20 bit left key . 20 bit right key (PC2)
sub compress($) {
    my $v   = shift;
    my $ret = "";
    for my $o (@compress_lut) {
        $ret .= substr( $v, $o, 1 );
    }
    return $ret;
}

# lut for expanding 20 bit data to make 30 bits (E)
my @extend_lut = ( 19, 4, 3, 8, 7, 12, 11, 16, 15, 0 );

# expand 20 bit right data to make 30 bits to xor with subkey
sub extend($) {
    my $v   = shift;
    my $ret = "";
    for my $o (@extend_lut) {
        $ret .= substr( $v, $o, 1 );
    }
    return $ret . $v;
}

# xor two equal length bitstrings
sub xors($$) {
    my ( $a, $b ) = @_;

    return "" . ( ( $a ^ $b ) | ( "\x30" x length($a) ) );
}

# lut for permutation after xor but before sboxes p1 (not in DES)
my @p1box_lut = (
    8,  9,  26, 27, 28, 29, 6,  7,  22, 23, 24, 25, 4,  5, 18, 19,
    20, 21, 2,  3,  14, 15, 16, 17, 0,  1,  10, 11, 12, 13,
);

# permute output of xor into 5x6 bit groups for the sboxes (not in DES)
sub p1box($) {
    my $v   = shift;
    my $ret = "";
    for my $o (@p1box_lut) {
        $ret .= substr( $v, $o, 1 );
    }
    return $ret;
}

# luts for the 5 sboxes
my @sbox_lut1 = (
    '1011', '1110', '0010', '0101', '0011', '1111', '0111', '1101',
    '0000', '1001', '0001', '1000', '1100', '0100', '1010', '0110',
    '0100', '1101', '1000', '1011', '0001', '1010', '0011', '0010',
    '1100', '0101', '0111', '0000', '1001', '1111', '0110', '1110',
    '0001', '0111', '1101', '0011', '1010', '0000', '1100', '1111',
    '0100', '1000', '0110', '1011', '0010', '1001', '1110', '0101',
    '1111', '0000', '0110', '0011', '0011', '1011', '1100', '1110',
    '0001', '0010', '1010', '0100', '0101', '1000', '1001', '0111',
);

my @sbox_lut2 = (
    '1011', '0000', '0010', '1100', '0111', '1101', '1111', '0110',
    '1001', '1110', '1010', '0011', '0001', '1000', '0100', '0101',
    '0110', '1000', '0010', '0000', '1101', '1011', '1001', '0101',
    '0001', '1100', '1010', '1111', '1110', '0111', '0011', '0100',
    '1100', '0111', '0011', '1011', '1010', '1110', '0000', '0110',
    '0001', '1000', '1111', '0010', '1101', '1001', '0101', '0100',
    '0001', '0111', '1001', '1101', '1011', '0000', '1111', '1110',
    '1000', '0100', '1100', '0011', '1010', '0101', '0010', '0110',
);

my @sbox_lut3 = (
    '1011', '1101', '1100', '0100', '0101', '0110', '1110', '1111',
    '0111', '0001', '1000', '1001', '0011', '0010', '1010', '0000',
    '0101', '1011', '0010', '1000', '1001', '0111', '0001', '1101',
    '0011', '1010', '0000', '1111', '0110', '0100', '1100', '1110',
    '0111', '1100', '0011', '1111', '1000', '0110', '1011', '0101',
    '0010', '1010', '1001', '1101', '0000', '0001', '0100', '1110',
    '1110', '0001', '0111', '0000', '1000', '0010', '1010', '1111',
    '0101', '1011', '1100', '1101', '0100', '0011', '1001', '0110',
);

my @sbox_lut4 = (
    '1010', '0011', '1111', '0111', '0110', '1110', '0000', '0100',
    '1001', '0001', '1101', '0101', '1000', '1100', '1011', '0010',
    '1100', '0001', '1101', '0101', '0100', '0011', '0000', '0111',
    '1001', '1111', '1110', '0010', '0110', '1011', '1010', '1000',
    '1011', '1100', '0100', '1111', '0101', '0110', '1110', '0011',
    '1000', '0010', '1001', '1101', '0000', '0111', '0001', '1010',
    '0100', '1011', '1111', '0111', '0001', '0101', '1010', '1001',
    '1101', '0110', '1100', '0011', '0000', '1000', '1110', '0010',
);

my @sbox_lut5 = (
    '1010', '0010', '0000', '1111', '0110', '0111', '1101', '1000',
    '0011', '1100', '1011', '0101', '1001', '0001', '0100', '1110',
    '0010', '1001', '0101', '1101', '1100', '1110', '1111', '1000',
    '0110', '0111', '1011', '0001', '0000', '1010', '0100', '0011',
    '1000', '0000', '1101', '1111', '0001', '1100', '0011', '0110',
    '1011', '0100', '1001', '0101', '1010', '0111', '0010', '1110',
    '0011', '1101', '0000', '1100', '1001', '0110', '1111', '1011',
    '0001', '1110', '1000', '1010', '0010', '0111', '0100', '0101',
);

# run the 30 bit output from the p1 permutation through
# the five s boxes to generate a 20 bit value

sub sboxen($) {
    my $v   = shift;
    my $ret = "";

    $ret .= $sbox_lut1[ be_bin_to_int( substr( $v, 0,  6 ) ) ];
    $ret .= $sbox_lut2[ be_bin_to_int( substr( $v, 6,  6 ) ) ];
    $ret .= $sbox_lut3[ be_bin_to_int( substr( $v, 12, 6 ) ) ];
    $ret .= $sbox_lut4[ be_bin_to_int( substr( $v, 18, 6 ) ) ];
    $ret .= $sbox_lut5[ be_bin_to_int( substr( $v, 24, 6 ) ) ];

    return $ret;
}

# lut for p2 permutation (P)

my @p2box_lut =
  ( 19, 3, 6, 8, 0, 15, 18, 13, 10, 16, 2, 1, 4, 9, 7, 11, 5, 12, 14, 17, );

# permute the output of the sbox stage to get the xorand

sub p2box($) {
    my $v   = shift;
    my $ret = "";

    for my $o (@p2box_lut) {
        $ret .= substr( $v, $o, 1 );
    }

    return $ret;
}

# decrypt the data using the not-quite-des algorythm

sub decrypt($$) {
    my ( $data, $key ) = @_;

    # print "k:" . b4($key) . "\n";

    my $kl = reverse( substr( $key, 0,  20 ) );
    my $kr = reverse( substr( $key, 20, 20 ) );

    $kl = ror($kl);
    $kr = ror($kr);

    my $dl = reverse( substr( $data, 0,  20 ) );
    my $dr = reverse( substr( $data, 20, 20 ) );

    # 16 rounds, DES numbers them backwards for decryption

    for ( my $r = 16 ; $r > 0 ; --$r ) {

        $kl = ror($kl);
        $kr = ror($kr);

        # rounds 8,7 and 3 ror the key twice.
        if ( ( $r == 8 ) || ( $r == 7 ) || ( $r == 3 ) ) {
            $kl = ror($kl);
            $kr = ror($kr);
        }

        # print "Round $r\n";
        # print "  kl:", b4($kl), "\n";
        # print "  kr:", b4($kr), "\n";

        my $k = compress( $kl . $kr );

        my $dre = extend($dr);

        # print "  k: ", b6($k), "\n";
        # print "  dl:",  b4($dl),  "\n";
        # print "  dr:",  b4($dr),  "\n";
        # print "  dre:", b6($dre), "\n";

        my $xor = xors( $dre, $k );

        # print "  xor:", b6($xor), "\n";

        my $p1 = p1box($xor);

        # print "  p1 :", b6($p1), "\n";

        my $s = sboxen($p1);

        # print "  sout:", b4($s), "\n";

        my $p2 = p2box($s);

        # print "  p2 :", b4($p2), "\n";

        my $ndl = $dr;
        my $ndr = xors( $dl, $p2 );

        # print "  ndl: ", b4($ndl),"\n";
        # print "  ndr: ", b4($ndr),"\n";

        $dl = $ndl;
        $dr = $ndr;
    }

    my $ck = substr( $dl, 0, 16 );

    # a valid weather packet always has this value

    return undef unless $ck == '0010010100000001';

    my $out = "";

    # one final gotcha the data are bit LE, but word BE.

    $out .= reverse( substr( $dr, 16, 4 ) );
    $out .= reverse( substr( $dr, 12, 4 ) );
    $out .= reverse( substr( $dr, 8,  4 ) );
    $out .= reverse( substr( $dr, 4,  4 ) );
    $out .= reverse( substr( $dr, 0,  4 ) );
    $out .= reverse( substr( $dl, 16, 4 ) );

    return $out;
}

######################################################### DECODE

sub time_to_string($$) {
    my ( $dt, $offset ) = @_;
    $dt = $dt->clone->add( days => $offset );
    return $dt->dmy('/');
}

my @day_lut = (
    ['--'],
    [ 'sunny', 'heatwave' ],
    ['slightly cloudy'],
    ['mostly cloudy'],

    ['overcast'],
    ['thunderstorms'],
    [ 'heavy rain', 'extreme rain' ],
    [ 'snow',       'heavy snow' ],

    [ 'fog',         'visbility below 50m' ],
    [ 'sleet',       'heavy sleet' ],
    [ 'rain shower', 'heavy rain showers' ],
    ['light rain'],

    [ 'snow showers', 'heavy snow showers' ],
    [ 'front storm',  'strong thunderstorms' ],
    ['high altitude fog'],
    [ '?sleet showers', '?heavy sleet showers' ],
);

my @night_lut = (
    ['--'],
    ['clear'],
    ['slightly cloudy'],
    ['mostly cloudy'],

    ['overcast'],
    ['thunderstorms'],
    [ 'heavy rain', 'extreme rain' ],
    [ 'snow',       'heavy snow' ],

    [ 'fog',         'visbility below 50m' ],
    [ 'sleet',       'heavy sleet' ],
    [ 'rain shower', 'heavy rain showers' ],
    ['light rain'],

    [ 'snow showers', 'heavy snow showers' ],
    [ 'front storm',  'strong thunderstorms' ],
    ['high altitude fog'],
    [ '?sleet showers', '?heavy sleet showers' ],
);

sub day($$) {
    my ( $v, $h ) = @_;

    my $ret = $day_lut[$v]->[$h];

    return $ret if defined $ret;

    $ret = "HEAVY(" . $day_lut[$v]->[0] . ")";

    return $ret;
}

sub night($$) {
    my ( $v, $h ) = @_;

    my $ret = $night_lut[$v]->[$h];

    return $ret if defined $ret;

    $ret = "HEAVY(" . $night_lut[$v]->[0] . ")";

    return $ret;
}

my @heavy_lut = (
    'no heavy weather',
    'heavy weather for 24hrs',
    'heavy weather in daytime',
    'heavy weather in nighttime',

    'storm for 24hrs',
    'storm in daytime',
    'storm in nighttime',
    'gusts in daytime',

    'gusts in nighttime',
    'freezing rain during morning',
    'freezing rain during afternoon',
    'freezing rain during night',

    'dust',
    'ozone',
    'radiation',
    'flood'
);

my @sun_lut = ( '0-2 hours', '2-4 hours', '5-6 hours', '7-8 hours' );
my @rain_lut = ( '0%', '15%', '30%', '45%', '60%', '75%', '90%', '100%' );

my @wind_dir_lut = (
    'N', 'NE', 'E', 'SE', 'S', 'SW', 'W', 'NW', '0001', '1001',
    '0101 Bise N/E',
    '1101 Mistral N',
    '0011', '1011', '0111', '1111'
);

my @wind_speed_lut =
  ( '0 Bft', '0-2 Bft', '3-4 Bft', '5-6 Bft', '7 Bft', '101', '011', '111' );

sub decode_bits($$) {
    my ( $bits, $odd ) = @_;

    my $day_val   = le_bin_to_int( substr( $bits, 0, 4 ) );
    my $night_val = le_bin_to_int( substr( $bits, 4, 4 ) );

    my $ret = {};

    my $heavy_day   = 0;
    my $heavy_night = 0;

    my $anomaly = int( substr( $bits, 15, 1 ) );

    $ret->{anomaly} = $anomaly;

    if ( $anomaly eq 0 ) {

        if ( $odd eq 0 ) {

            my $heavy_val = le_bin_to_int( substr( $bits, 8, 4 ) );

            $ret->{heavy} = $heavy_lut[$heavy_val];

            if ( $heavy_val eq 1 ) {
                $heavy_day   = 1;
                $heavy_night = 1;
            }
            elsif ( $heavy_val eq 2 ) {
                $heavy_day = 1;
            }
            elsif ( $heavy_val eq 3 ) {
                $heavy_night = 1;
            }

        }
        else {

            my $wind_dir_val   = le_bin_to_int( substr( $bits, 8,  4 ) );
            my $wind_speed_val = le_bin_to_int( substr( $bits, 12, 3 ) );

            $ret->{wind_speed} = $wind_speed_lut[$wind_speed_val];
            $ret->{wind_dir}   = $wind_dir_lut[$wind_dir_val];

        }

    }
    else {
        my $morning_val = le_bin_to_int( substr( $bits, 8,  2 ) );
        my $sun_val     = le_bin_to_int( substr( $bits, 10, 2 ) );

        $ret->{morning} = "morning " . $morning_val;
        $ret->{sun}     = $sun_lut[$sun_val];
    }

    if ( $odd eq 0 ) {
        my $rain_val = le_bin_to_int( substr( $bits, 12, 3 ) );
        $ret->{rain} = $rain_lut[$rain_val];
    }

    $ret->{day} = day( $day_val, $heavy_day );
    $ret->{night} = night( $night_val, $heavy_night );

    my $temp_val = le_bin_to_int( substr( $bits, 16, 6 ) );
    if ( $temp_val == 0 ) {
        $ret->{temp} = "< -21⁰C";
    }
    elsif ( $temp_val == 0x3f ) {
        $ret->{temp} = "> 41⁰C";
    }
    else {
        $ret->{temp} = ( -22 + $temp_val ) . "⁰C";
    }

    return $ret;
}

sub print_prognosis($) {
    my $w = shift;
    print "  day        = " . $w->{day} . "\n";
    print "  night      = " . $w->{night} . "\n";
    print "  temp       = " . $w->{temp} . "\n";
}

sub print_rest($) {
    my $w = shift;
    print "  anomaly    = " . $w->{anomaly} . "\n";
    print "  heavy      = " . $w->{heavy} . "\n" if exists $w->{heavy};
    print "  rain       = " . $w->{rain} . "\n"  if exists $w->{rain};
    print "  wind speed = " . $w->{wind_speed} . "\n"
      if exists $w->{wind_speed};
    print "  wind dir   = " . $w->{wind_dir} . "\n" if exists $w->{wind_dir};
    print "  morning    = " . $w->{morning} . "\n"  if exists $w->{morning};
    print "  sun        = " . $w->{sun} . "\n"      if exists $w->{sun};
}

my @location_lut = (
    "FR Bordeaux",
    "FR la Rochelle",
    "FR Paris",
    "FR Brest",
    "FR Clermont-Ferrand",
    "FR Beziers",
    "BE Bruxelles",
    "FR Dijon",
    "FR Marseille",
    "FR Lyon",
    "FR Grenoble",
    "CH La Chaux de Fonds",
    "DE Frankfurt am Main",
    "DE Trier",
    "DE Duisburg",
    "GB Swansea",
    "GB Manchester",
    "FR le Havre",
    "GB London",
    "DE Bremerhaven",
    "DK Herning",
    "DK Arhus",
    "DE Hannover",
    "DK Copenhagen",
    "DE Rostock",
    "DE Ingolstadt",
    "DE Munich",
    "IT Bolzano",
    "DE Nuernberg",
    "DE Leipzig",
    "DE Erfurt",
    "CH Lausanne",
    "CH Zurich",
    "CH Adelboden",
    "CH Sion",
    "CH Glarus",
    "CH Davos",
    "DE Kassel",
    "CH Locarno",
    "IT Sestriere",
    "IT Milan",
    "IT Rome",
    "NL Amsterdam",
    "IT Genova",
    "IT Venezia",
    "DE Strasbourg",
    "AT Klagenfurt",
    "AT Innsbruck",
    "AT Salzburg",
    "AT/SK Vienna/Bratislava",
    "CZ Prague",
    "CZ Decin",
    "DE Berlin",
    "SE Gothenburg",
    "SE Stockholm",
    "SE Kalmar",
    "SE Joenkoeping",
    "DE Donaueschingen",
    "NO Oslo",
    "DE Stuttgart",
    "IT Naples",
    "IT Ancona",
    "IT Bari",
    "HU Budapest",
    "ES Madrid",
    "ES Bilbao",
    "IT Palermo",
    "ES Palma de Mallorca",
    "ES Valencia",
    "ES Barcelona",
    "AD Andorra",
    "ES Sevilla",
    "PT Lissabon",
    "IT Sassari",
    "ES Gijon",
    "IE Galway",
    "IE Dublin",
    "GB Glasgow",
    "NO Stavanger",
    "NO Trondheim",
    "SE Sundsvall",
    "PL Gdansk",
    "PL Warsaw",
    "PL Krakow",
    "SE Umea",
    "SE Oestersund",
    "CH Samedan",
    "HR Zagreb",
    "CH Zermatt",
    "HR Split"
);

my @cycle_4_lut = (
    "Highs day 1",
    "Lows day 1",
    "Highs day 2",
    "Lows day 2",
    "Highs day 3",
    "Lows day 3",
    "Highs day 4",
    "Wind day 4",
);

my @cycle_2_lut = ( "Forecast day 1", "Forecast day 2", );

sub decode($$$) {
    my ( $dt, $mins_since_midnight, $bits ) = @_;

    $mins_since_midnight = int( ( $mins_since_midnight - 1 ) / 3 );

    my $cycle  = int( $mins_since_midnight / 60 );
    my $region = $mins_since_midnight % 60;
    my $day    = int( $cycle / 2 );

    my $w = decode_bits( $bits, $cycle % 2 );

    if ( $cycle < 7 ) {
        print "region: ", $location_lut[$region], ", cycle: ",
          $cycle_4_lut[$cycle], ", when: ", time_to_string( $dt, $day ), "\n";
        print_prognosis($w);
        print_rest($w);
    }
    else {

        print "region: ", $location_lut[$region], ", cycle: ",
          $cycle_4_lut[$cycle], ", when: ", time_to_string( $dt, $day ), "\n";
        print_rest($w);

        my $sub_region = $region % 30;
        my $sub_cycle  = int( $region / 30 );
        my $sub_day    = 1 + $sub_cycle;

        $sub_region += 60;

        print "region: ", $location_lut[$sub_region], ", cycle: ",
          $cycle_2_lut[$sub_cycle], ", when: ",
          time_to_string( $dt, $sub_day ), "\n";
        print_prognosis($w);

    }

}

######################################################### MAIN

my $d = [];
my $n = 0;

open DATA, "<", "raw";

while (<DATA>) {
    chomp;
    $d->[$n] = $_;
    $n++;
}

for ( my $i = 0 ; $i < ( $n - 2 ) ; ++$i ) {
    my ( $data_bits, $key_bits,$cest ) =
      extract( $d->[$i], $d->[ $i + 1 ], $d->[ $i + 2 ] );

    my $min   = le_bcd_to_int( substr( $key_bits, 0,  8 ) );
    my $hour  = le_bcd_to_int( substr( $key_bits, 8,  8 ) );
    my $day   = le_bcd_to_int( substr( $key_bits, 16, 8 ) );
    my $month = le_bcd_to_int( substr( $key_bits, 24, 5 ) );
    my $year  = le_bcd_to_int( substr( $key_bits, 32, 8 ) );
    $year += 2000;

    my $mins_since_midnight = $min + ( $hour * 60 );

    next unless ( $mins_since_midnight % 3 ) == 2;

    $mins_since_midnight+=60 if $cest == "0";

    $mins_since_midnight-=1440 while $mins_since_midnight>1440;



    my $plain = decrypt( $data_bits, $key_bits );

    next unless defined $plain;

    my $dt = DateTime->new(
        year  => $year,
        month => $month,
        day   => $day
    );

    printf "%04d-%02d-%02d %02d:%02d line %8d bits: %s\n", $year, $month, $day,
      $hour, $min, $i, b4($plain);

    decode( $dt, $mins_since_midnight, $plain );
}