Dotaz: Jak zjistit počet výskytů 16b kombinací v Python3 ?

[Zm]

Ahoj, chci číst binární soubor po dvou bajtech a u každé z 65536 kombinací zaznamenat počet výskytů. Na konci programu vypsat, která z kombinací se vyskytuje méně než nějaká hranice (např. 100x).

V C bych si udělal uint64_t pole[65536] a po načtení 2B inkrementoval hodnotu. Jak to udělat v Python 3 aby to bylo správně pythonistické?

---

---

Odpovědi

[Zm]

V C to mám napsané takto

#include string.h>
#include stdio.h>

#define COMB_COUNT 65536

int main(void)
{
  const char *filename = "file.bin";

  unsigned long long int treshold = 5;
  unsigned long long int bin_arr[COMB_COUNT];
  unsigned twob;
  unsigned int i;
  FILE *fd;

  if (!(fd = fopen(filename, "rb"))) {
    return -1;
  }

  memset(bin_arr, '\0', COMB_COUNT);
  while ((i = fread(&twob, 2, 1, fd)) > 0) {
    bin_arr[twob] += 1;
  }

  for (i = 0; i < COMB_COUNT; i++) {
    if (bin_arr[i] < treshold) {
      printf("%04x, %llu times\n", i, bin_arr[i]);
    }
  }

  fclose(fd);

  return 0;
}

[Gréta]

muj pokus o co nejpodobnější přepsání do pythonu :O :D :D ;D

#!/usr/bin/env python3

POCET_KOMBINACI = 65536
TRESHOLD = 100

try:
    fd = open('soubor.bin', 'rb')
except IOError:
    print('nejde otevrit soubor binarni jeden :O :\'(') # :D :D
    exit(-1)

with fd:
    # uděláme si list vo velikosti počtu kombinací
    bin_arr = [0] * POCET_KOMBINACI
    
    # zkusíme přečíst dva byty
    byty = fd.read(2)
    
    # dokud máme dva byty přečtený
    while len(byty) == 2:
    
        # převedem na index v tom našem listu 'bin_arr'
        # v metodě from_bites musíme říct jakej to je endian jakože jestli big nebo little
        index = int.from_bytes(byty, 'little') 
        bin_arr[index] += 1
        
        # načtem dvojici bytů pro další iteraci
        byty = fd.read(2)
        
    # zavřem soubor
    fd.close()
        
    # vypišem ty věčí hodnoty než treshold
    for i in range(POCET_KOMBINACI):
        if bin_arr[i] > TRESHOLD:
            print(f'{i.to_bytes(2, "little")}, {bin_arr[i]} times')
            
    # v pythonu mužem je asi obyklejší iterovat foreachem
    # takže by to asi jako lepší bylo
    # for hodnota in bin_arr:
    #    if hodnota > TRESHOLD:
    #        ....
            
    exit(0) # zbytečný volat nakonec skriptu v pythonu si myslim ale nevim jistě :O :O

[MadCatX]

import sys
import struct

log = {}
with open(sys.argv[1], 'rb') as fh:
    raw = fh.read(2)
    while raw is not None and len(raw) == 2:
        word = struct.unpack('<H', raw)[0]
        if word in log:
            log[word] += 1
        else:
            log[word] = 1
        raw = fh.read(2)

for k in sorted(log.keys()):
    v = log[k]
    if v < 100:
        print('[{}]: {}'.format(k, v))

[Jendа]

        if word in log:
            log[word] += 1
        else:
            log[word] = 1

log = defaultdict(int)

[Jendа]

1) int pole[100] = {0}; // immediate inicializace

2)

c.c: In function ‘main’:
c.c:20:3: warning: ‘memset’ used with length equal to number of elements without multiplication by element size [-Wmemset-elt-size]
   20 |   memset(bin_arr, '\0', COMB_COUNT);
      |   ^~~~~~

Doporučuju používat nějaký hustý kompilátor (třeba rozumně nové gcc) se zapnutými všemi warningy (tohle bylo -Wall -Wextra -Wduplicated-cond -Wduplicated-branches -Wlogical-op -Wrestrict -Wnull-dereference -Wjump-misses-init -Wdouble-promotion -Wshadow -Wformat=2).

[Zm]

Nebo tam dát

memset(bin_arr, '\0', sizeof(bin_arr));

[z_sk]

Lepsie by bolo, keby twob bol typ uint16_t, vid..

[Gréta]

atak já ti to teda zkusim napsat ale jako fakt nevim jestli je to jakože bude fakt správnej pythonistickej zápis :O ;D

# nóó začala bych jakoby tak že bych si udělala slovník kterej bude mit
# jako keys/klíče ty jednotlivý kombinace bytů a jako values/hodnoty
# počty nalezenejch výskytů 

# takže uděláme slovník
slovnik = {}

# a do něj nastrkáme ty klíče pro všecky ty konbinace dvojic bytů
for i in range(256):
    for j in range(256):
        
        # do slovníku přidáme novej klíč a hodnotu normálně jako slovnik[klic]=hodnota
        # tamto bytes([nejakej_int]) předělá int na bytes
        # noa byty mužemenormálně strkat zasebe tim plusem
        
        # tady ta podmínka if i > j je tady proto abysme ždycky ten věčí byte strčili 
        # nazačátek tý naší kombinace protože 'ab' je pronás to samý co 'ba'
        # neni pro todleto náhodou nějakej víc lepšejší matematickej způsob????? :O :O
        if i > j:
            slovnik[ bytes([i]) + bytes([j])] = 0
        else:
            slovnik[ bytes([j]) + bytes([i])] = 0

        # poznámka:
        # nóó jak si jako teďko líp čtu tu tvou votázku :D tak jestli todleto nepotřebuješ a chceš jakože fakt kombinace jenom po sobě jdoucích bytů
        # tak ty ify dej do pryč a v ckylu postupuj jenom po jednom bytu a pamatuj si ten zminulý iterace a znich si dělej ty kombinace :O ;D

# čteme soubor binárně
with open('soubor.bin', 'rb') as f:
    
    # přečtem dvojici bytů takle
    prvni_byte = f.read(1)
    druhy_byte = f.read(1)
    # pokud byte neni read vrátí b''
    
    # TODO ošetřit situaci kdy máš lichej počet bytů :O :O :O :O
    # b'' se vyhodnotí jako false všecko vostatní true
    # nóó takže dokavaď máme dva byty
    while prvni_byte and druhy_byte:
        
        # kouknem jestli je první věčí než druhej jestli ne tak prohodíme
        # znova pokavaď todleto nepotřebuješ dej do pryč
        if druhy_byte > prvni_byte:
            tmp = prvni_byte
            prvni_byte = druhy_byte
            druhy_byte = tmp
          
        #sestavíme tu kombinaci
        kombinace = prvni_byte + druhy_byte
        # a inkrementujem hodnotu ve slovníku
        slovnik[ kombinace ] += 1
        
        # načteme byty pro další iteraci
        prvni_byte = f.read(1)
        druhy_byte = f.read(1)
        
        
# pro vokrasu seřadíme podle počtu výskytů těch kombinací
serazeno = sorted(slovnik.items(), key=lambda x: x[1])

# a vypišem všecky kde je víc jak těch tvejch 100 kombinací :O ;D
for klic, hodnota in serazeno:
    if hodnota > 100:
        print(f"kombinace: {klic} pocet: {hodnota}")

print('hotovo!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!')

[Zm]

Díky!

Trochu jsem doufal, že to bude v Pythonu efektivnější ve smysl počtu řádku, ale je to cca 20 C a 24 py, byť to py navíc seřadí.

Ale moc se mi nezdá to prohazování bajtů podle velikosti. Pokud to beru tak, že mám všechny kombinace 16 bitů a chci vědět kolikrát je která zastoupená, tak pro mě 'ab' není to samý co 'ba'.

[Gréta]

nóó když neni tak neni :D ;D btw v komentu v poznámce pišu jestli to jako nechceš tak abys jakoby umazal ty ify a prohazování ;D

pokusila sem se ještě přepsat ten tvuj cčkovej zdrojáček do pythonu tak koukni jestli to je jakoby to co ty chceš :O ;D

[Zm]

Napadlo mě jestli spíš nemám řešit endianitu. Ale vzhledem k tomu, že se to změnou endianu přesune jen na jiné místo v tom poli, je to asi zbytečné.

Každopádně velké dík za pomoc.

[Bherzet]

Ale vzhledem k tomu, že se to změnou endianu přesune jen na jiné místo v tom poli, je to asi zbytečné.

Ano, ke kolizím docházet nebude, ale vypíše to jiný výsledek (word 0xaabb na vstupu to bude prezentovat jako 0xbbaa na výstupu).

[Gréta]

podobně jako v madcataxovým zdrojáčku de tady taky nahradit ten slovnik={} a jeho 'alokace' tim jendovým defaultdictem slovnik=defaultdict(int)

[Jendа]

Tak schválně jak si ta řešení stojí z hlediska výkonu. Hledáš věci, které se vyskytují alespoň 100x, takže to asi nebude úplně malé. Osobně tipuju, že zde uvedené Pythoní řešení budou 100x pomalejší než naivní C, a že použitím nedávno zde probíraného fread_unlocked to C zrychlím ještě 2x.

Tak schválně.

Pro tento účel jsem pořídil 19MB náhodný soubor, kde se vyskytovala jedna kombinace 98x, jinak všechny víckrát.

Naivní C s ASANem: 0.25s

Naivní C -O3: 0.21s

Naivní C -O3 unlocked: 0.17s

Gréta 1 (otočil jsem jí nerovnost): 2.9s

MadCatX: 4.4s

MadCatX s defaultdict: 4s

Závěr: není to tak špatné jak jsem si myslel, ale furt je to pomalejší cca. 20x.

Osobně bych to dal jako modul a natáhl ho přes CFFI. Chci o tom napsat blogpost.

[MadCatX]

mmapované C vs Python:

C:

#include <stdio.h>
#include <stdlib.h>
#include <sys/mman.h>
#include <sys/types.h>
#include <sys/stat.h>
#include <assert.h>
#include <fcntl.h>
#include <unistd.h>
#include <stdint.h>
#include <string.h>

#define NUM_WORDS 65536
#define RESULTS_SIZE (sizeof(uint_fast64_t) * NUM_WORDS)

int main(int argc, char **argv)
{
        int fd;
        void *data;
        size_t data_length;
        unsigned short *pos;
        unsigned short *end;
        uint_fast64_t *results;

        if (argc < 1)
                return EXIT_FAILURE;

        results = malloc(RESULTS_SIZE);
        if (results == NULL)
                return EXIT_FAILURE;

        memset(results, 0, RESULTS_SIZE);

        fd = open(argv[1], O_RDONLY);
        if (fd < 0)
                return EXIT_FAILURE;

        data_length = lseek(fd, 0, SEEK_END);
        data_length -= data_length % 2;
        data = mmap(NULL, data_length, PROT_READ, MAP_PRIVATE, fd, 0);
        if (data == MAP_FAILED)
                return EXIT_FAILURE;

        pos = (unsigned short *)data;
        end = (unsigned short *)(data + data_length);
        for (; pos < end; pos++)
                results[*pos]++;

        for (size_t idx = 0; idx < NUM_WORDS; idx++) {
                uint_fast64_t value = results[idx];
                if (value < 1500)
                        printf("[%zu]: %lu\n", idx, value);
        }

        return EXIT_SUCCESS;
}

Totéž v Pythonu:

import mmap
import struct
import sys


log = [0] * 65536
with open(sys.argv[1], 'rb') as fh:
    fh.seek(0, 2)
    data_length = fh.tell()
    data_length -= data_length % 2
    data = mmap.mmap(fh.fileno(), data_length, mmap.MAP_PRIVATE, mmap.PROT_READ)

    unpack_word = struct.Struct('<H').unpack_from
    idx = 0
    while idx < data_length:
        word = unpack_word(data, idx)[0]
        log[word] += 1
        idx += 2

for k,v in enumerate(log):
    v = log[k]
    if v < 200:
        print('[{}]: {}'.format(k, v))

Časy pro 20, resp. 200 MB soubory z /dev/urandom:

C: 0,22; 0,35 s

Python: 4,8; 49,0 s

U Cčka je dost znát ten fixní overhead, zatímco Python visí v té čtecí smyčce. No co, aspoň je vidět, že na různé problémy je vhodné vzít si vždy správný šroubovák.

[Jendа]

Tohle je ještě 3x rychlejší než to tvoje:

import numpy as np

a = np.fromfile("r", dtype=np.uint16)
bins = np.arange(0, 65535+1, 1)
hist, _ = np.histogram(a, bins)

ma_data = np.ma.masked_greater(hist, 100)

print(ma_data.compressed())

[MadCatX]

Hustý, něco takového jsem hledal, když jsem vymýšlel, jak na to napasovat Numbu.

import sys
import numba
import numpy


@numba.jit(nopython=True)
def walk(buf):
    log = [0] * 65536
    for v in buf:
        log[v] += 1
    return log


with open(sys.argv[1], 'rb') as fh:
    data = numpy.fromfile(fh, numpy.uint16)
    log = walk(data)

for k, v in enumerate(log):
    v = log[k]
    if v < 150:
        print('[{}]: {}'.format(k, v))

Dává do docela zajímavé výsledky:

20 MB, práh 150:
NumPy: 1,15 s
NumPy+Numba: 1,00 s

200 MB, práh 1500:
NumPy: 10,2 s
NumPy+Numba: 1,25 s

Aneb nejrychlejší způsob, jak zrychlit Python je zbavit se co nejvíce nativně pythoního kódu :)

[Jindřich Makovička]

Případně jestli rádi Pandas...

#!/usr/bin/python3

import numpy as np
import pandas as pd

pd.set_option("display.max_rows", None)

dt = np.dtype([('value', 'u2')])
data = np.fromfile('data.bin', dtype=dt)
df = pd.DataFrame(data)

df2 = df.groupby(by='value')['value'].count()

print(df2[df2 < 10])

[Bherzet]

mmap v Pythonu jsem zkoušel včera a vyšlo mi to na cca 150MB soubor skoro 2× pomalejší než prostý read po bajtech.

[Zm]

Ten binární soubor má okolo 100 MB.

[Jendа]

A vadí když to poběží 20 sekund? (ty moje testy byly na notebooku s i5 8. generace, takže poměrně moderní, ale notebook)

Nechal bych to v C a přilepil přes CFFI.

[Zm]

Pouštím to ručně, takže bych těch 20s přežil. Ale čekat do 1s v C je příjemnější.

[MadCatX]

I v pythonu se dostaneš k časům okolo 1 s, když ty nejprotěžovanější cykly necháš JITnout Numbou. Zrovna pro tvůj případ se to zdá jako šikovné řešení.

[Gréta]

nóó jestli to jako neni děsně tajný proč to jakoby potřebuješ mit napsaný v tom pythonu misto v cčku?? :O :O

[Zm]

Není na tom nic tajného. Pomocné prográmky si raději dělám v Pythonu, protože se v něm většínou udělají snadněji (tento případ je, zdá se, výjimka potvrzující pravidlo), nemusí se kompilovat, jsou blbuvzdornější (viz. to pomotané memset) a taky proto abych se Py. trochu naučil.

[Gréta]

nóóóó ještě de udělat hulvátřešení jakože zkompilovat si cčkovej kód do binárky atu pouštět z pythonu a výstup si v pouštěcím skriptu parsovat :O ;D čistý ani pythonistický řešení to asi jako neni ale funguje to :D ;D tamto cffi asi jako bude vo něco rychleší si myslim

[Gréta]

Chci o tom napsat blogpost.

hele už to napsal ;D

[Andrej]

import sys
if len(sys.argv) != 2:
  print('Hovno, vole.', file=sys.stderr)
  sys.exit(3)
pairs = {}
with open(sys.argv[1], 'rb') as fuck:
  while shit := fuck.read(2):
    if len(shit) == 1:
      print('Zasraný soubor.', file=sys.stderr)
      sys.exit(5)
    pairs.setdefault(shit, [0])[0] += 1
counts = {}
for pair, count in pairs.items():
  counts.setdefault(count[0], []).append(pair)
for key in sorted(counts.keys(), reverse=True):
  for b in sorted(counts[key]):
    print(f'{key}: [{b[0]:02x}, {b[1]:02x}]')

[BFU]

Kdyz ono to v shellu jde jako oneliner

$ xxd -c2 -ps file | sort | uniq -c | sort

[Gréta]

ti tam jakoby eště chybí ten threshold ;D

[BFU]

No tak zrovna pro min nez 100 vyskytu treba grep "^ \+[0-9][0-9]\? ", ackoli by me zajimalo, jestli je neco lepsiho nez to drtit regexama.

[Gréta]

možná jakoby ňák rozpůlit výstupy toho posledního sortu mezerou/bílejma znakama a první půlky předělat na int a koukat jestli sou věčí než treshold a pak to vypysovat ale toby jako byl mocmocmoc vošlivej oneliner :O :D :D ;D

[Charon]

Jo, tak ten hexdump je dost dobrý trik!!! To mne inspirovalo k úpravám v jiných skriptech, kde to "řeším" přes gawk. Takže moc díky za tip. Na oplátku dodám "improvement" navrhovaného one-lineru.

xxd -c2 -ps file | LANG=C sort --parallel=8 | uniq -c | sort -k 2 -k 1

Pro můj testovací 50MB file to na mém MacBook běží v parallel verzi cca 37sekund. Původně mě to běželo 3 minuty, což mě nepřijde jako "promptní" one-liner. A proto ta motivace k inovaci na rychlost a úsporu času. Může být zajímavé v sort použít přepínač -S 8G, kdy si příkaz sort alokuje více RAM paměti na obcování s daty při jejich třídění.

Řešení, jak vypsat výskyty menší než trashold? (např. méně než 100), nechám dalším ke studiu. Byť ten návrh s grep mě přijde OK.

Bye. Charon

[Gréta]

ešení, jak vypsat výskyty menší než trashold? (např. méně než 100), nechám dalším ke studiu. Byť ten n

xxd -c2 -ps file | LANG=C sort --parallel=8 | uniq -c | sort -k 2 -k 1 | awk -v thr="100" '$1 > thr'

[Jendа]

Může být zajímavé v sort použít přepínač -S 8G

A přitom ta úloha má zjevné lineární řešení s konstantní pamětí a sort třídí v nlogn čase a lineárně paměti

[Gréta]

nj ale když uděláš něco takovýho

xxd -c2 -ps file | awk '{a[$0]++}END{for(i in a){if(a[i]>100)print a[i], i}}' | sort

tak už to jakoby neni žádnej oneliner hezkej :/ :/