はじめに

JPEG File を DecodeすることによってJPEGを理解します。
JPEG は ITU-T T.81 で規格化されています。

Image component

画像全体をimageといいます。imageは画素で構成されます。
水平、垂直方向の画素数をそれぞれ$X,Y$と表記します。
画素はY'(輝度)とCb Cr(色差)で構成されます。Y',Cb,Crをimage componentといいます。

視覚はY'に比べてCb/Crの変化を感じにくいため、色差のsample数を減らすことができます。
Y'の水平、垂直方向のsample数はimageの画素数$X,Y$と同じです。
Y'のsample数に対するCb/Crのsample数をsampling factorで指定します。

Sampling factor

sampling factorを$H_i, V_i$と表記します。i=1がY',i=2がCb,i=3がCrです。
componentの水平、垂直方向のsample数($x_i, y_i$)は次の式で求めます。

x_i = \lceil X \frac{H_i}{H_{max}} \rceil \\
y_i = \lceil Y \frac{V_i}{V_{max}} \rceil \\

式1とします。$\lceil \rceil $ は整数への切り上げを表します。
$H_{max}, V_{max}$ はそれぞれ $H_i, V_i$ の最大値です。

具体例(例1)で説明します。
$H_i,V_i$が次のように与えられているとします。

Name	$H_i$	$V_i$
Y'(i=1)	2	2
Cb(i=2)	1	1
Cr(i=3)	1	1
$H_{max},V_{max}$	2	2

$X,Y$が次のように与えられているとします。

Name	Value
X	30
Y	31

式1から$x_i,y_i$は次の値になります。

Name	$x_i$	$y_i$
Y'(i=1)	$\lceil 30*\frac{2}{2} \rceil = 30$	$\lceil 31*\frac{2}{2} \rceil = 31$
Cb(i=2)	$\lceil 30*\frac{1}{2} \rceil = 15$	$\lceil 31*\frac{1}{2} \rceil = 16$
Cr(i=3)	$\lceil 30*\frac{1}{2} \rceil = 15$	$\lceil 31*\frac{1}{2} \rceil = 16$

Data unit

8x8 sample を data unitといいます。
Encoder / Decoder は data unit の単位でDCT(Discrete Cosine Transform)を行います。

Zig-zag order

EncoderはFDCTの後、sample 単位で Entropy coding をします。
そのため data unitに含まれる一つ一つの sample を順番に並べる必要があります。
この順番をZig-zag orderといいます。次に示す番号順です。0をDC, 1-63をACといいます。

Minimum coded unit (MCU)

Zig-zag orderはsampleの並び順でしたが、data unitについても並び順を考える必要があります。
component毎に$H_i $ x $V_i$の小さな矩形を考えます。
矩形に含まれる data unit を左上から右下に並べます。
それをcomponentの順番に並べてminimum coded unit(MCU)とします。

具体例(例2)を示します。sampling factor H1=2,V1=2 H2=1,V2=1 H3=1,V3=1 とします。
Y'について2x2の矩形を考え、左上から右下に向かって4つ並べます。
Cb,Crについて、1x1の矩形を考え、Cb data unitを1つ、その後ろにCrのdata unitを1つ並べます。
1つ目の MCU ができました。image全体ができるまでこれを繰り返します。

Encode / Decode は MCU 単位で行われます。
したがって Encode / Decodeされるcomponentの水平、垂直方向のsample数($x'_i, y'_i$)は次の式で与えられます。

x'_i = \lceil x_i, H_i*8 \rceil \\
y'_i = \lceil y_i, V_i*8 \rceil \\

式2とします。$\lceil a,b \rceil $ はaをbの倍数へ切り上げることを表します。

具体例として例1の数値の場合、$x'_i,y'_i$は次のようになります。

Name	$x'_i$	$y'_i$
Y'(i=1)	$\lceil 30, 16 \rceil = 32$	$\lceil 31, 16 \rceil = 32$
Cb(i=2)	$\lceil 15, 8 \rceil = 16$	$\lceil 16, 8 \rceil = 16$
Cr(i=3)	$\lceil 15, 8 \rceil = 16$	$\lceil 16, 8 \rceil = 16$

Encoder は右端、下端のsample が不足している場合、これを補います。

したがって Decoder は ECS の sample数が $x'_i,y'_i$ であると想定してDecodeできます。

Decode後に $x_i,y_i$ を有効sampleとして出力します。

example

以降の説明では次のJPEG Fileを具体例として使用します。earth.jpg とします。

earth.jpg の Hex dump は次の通りです。

$ od -vtx1 -Ax earth.jpg
000000 ff d8 ff e0 00 10 4a 46 49 46 00 01 01 01 00 48
000010 00 48 00 00 ff e1 00 22 45 78 69 66 00 00 4d 4d
000020 00 2a 00 00 00 08 00 01 01 12 00 03 00 00 00 01
000030 00 01 00 00 00 00 00 00 ff db 00 43 00 02 01 01
000040 02 01 01 02 02 02 02 02 02 02 02 03 05 03 03 03
000050 03 03 06 04 04 03 05 07 06 07 07 07 06 07 07 08
000060 09 0b 09 08 08 0a 08 07 07 0a 0d 0a 0a 0b 0c 0c
000070 0c 0c 07 09 0e 0f 0d 0c 0e 0b 0c 0c 0c ff db 00
000080 43 01 02 02 02 03 03 03 06 03 03 06 0c 08 07 08
000090 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c
0000a0 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c
0000b0 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c
0000c0 0c 0c ff c0 00 11 08 00 1f 00 1e 03 01 22 00 02
0000d0 11 01 03 11 01 ff c4 00 17 00 01 01 01 01 00 00
0000e0 00 00 00 00 00 00 00 00 00 00 09 08 0a 05 ff c4
0000f0 00 2f 10 00 02 01 03 03 03 03 01 07 05 00 00 00
000100 00 00 00 01 02 03 04 05 11 06 07 12 08 09 21 00
000110 13 41 0a 14 15 22 23 31 32 42 51 52 61 71 91 ff
000120 c4 00 17 01 00 03 01 00 00 00 00 00 00 00 00 00
000130 00 00 00 00 00 02 04 03 ff c4 00 1a 11 00 02 03
000140 01 01 00 00 00 00 00 00 00 00 00 00 00 00 01 02
000150 11 12 03 21 ff da 00 0c 03 01 00 02 11 03 11 00
000160 3f 00 83 7b 17 f6 2e bc 77 6c d7 57 5b ee a2 bb
000170 5c 34 8e d0 e9 1a a4 a3 ba dd 28 a2 0d 59 75 ab
000180 2a b2 1a 0a 32 e1 a3 59 16 26 57 92 46 57 11 09
000190 22 ca 39 90 0f 5a 31 d9 6e cc 9d 2e 74 e7 a3 29
0001a0 6c fa 67 61 76 ce 64 a5 ce 2b 6f 96 18 2f 97 19
0001b0 4b 7e e2 f5 55 8b 2c c7 3f da 18 28 fe 2a 07 8f
0001c0 5c ae c9 fb 35 66 e9 bf b5 86 c3 58 6c f0 ad 38
0001d0 ae d1 f6 fd 43 5b cb 1e e4 d5 b7 38 12 e1 50 ce
0001e0 7e 71 25 49 40 49 f0 91 a0 18 0a 31 58 d7 5e 92
0001f0 3a 4f ea 71 81 fe b1 e8 00 ef eb db b0 3f 4c 5d
000200 59 68 da aa 7a 7d b6 b0 ed 86 a1 8a 32 94 3a 8b
000210 44 50 c5 68 92 89 89 e5 ee 4b 4b 0a ad 2d 4a 02
000220 32 e2 48 fd c2 9c 82 4b 11 3c 86 67 fa ef e8 a3
000230 55 76 fc ea 52 fd b6 9a c3 ec b3 d6 da 5c 4d 45
000240 72 a6 52 68 ef 14 6f 9f 6a aa 03 f2 8d c5 81 1f
000250 aa 3a 48 87 f1 29 f5 b1 fd d3 d5 ad 6f a5 66 a7
000260 92 45 9b 89 5e 28 70 d8 23 07 18 ff 00 04 fe bf
000270 1f f4 03 ff 00 54 ee 8f b4 eb 9d af db 8d c2 a7
000280 89 20 bc 5a 2f 12 d8 6a 88 50 0c f1 d5 47 34 e9
000290 83 fa 94 43 46 c4 0f 01 5e 69 70 30 70 29 8f 0b
0002a0 e7 b1 75 ed 09 df 62 0e ae 6d 3d 50 f6 ad da 1b
0002b0 8d ba 68 3e f2 d0 d6 68 b4 2d f6 8e 30 5e 4a 4a
0002c0 9b 64 11 41 19 60 d8 c9 96 91 69 e7 c8 f1 f9 d8
0002d0 c8 20 81 66 d2 de 23 b5 d2 cd f7 94 94 f2 47 4a
0002e0 39 bc 43 1c 8a 1c 85 60 49 5e 20 90 47 93 8f 58
0002f0 ec ed 9d dc e3 5f f6 bf de d3 aa b4 6c b0 dc ac
000300 f7 73 0d 3e a2 d3 b5 b2 3a d1 5f 61 8c 92 80 95
000310 f3 14 f1 f3 73 14 ea 09 4f 71 d4 87 8e 49 62 91
000320 be d9 2f aa 73 a6 bd 7b a0 e9 64 bf 5c 2f 5b 73
000330 71 a7 8a 36 6b 6d ce c9 53 58 90 ca 14 82 21 96
000340 92 39 c3 aa e4 85 67 11 92 3c f0 53 e0 63 1a 68
000350 62 cf de fd cf b3 d4 57 dd 23 9a be 4b 0b 34 88
000360 96 f9 a5 e4 ca 64 e5 1a b0 05 55 b3 8e 4e 4f c0
000370 09 e4 8c f8 0b 7e a4 5d fc a0 ae a9 d3 bb 69 6b
000380 ad 82 e3 2c 97 69 75 2d 61 88 f2 f6 51 12 4a 6a
000390 71 9f 9c b4 95 43 3f 22 35 6f 21 81 34 17 5c bf
0003a0 52 86 d3 d1 d9 ab e1 da 7a 1a fd 6b a9 2a 62 f6
0003b0 e9 ab 24 a3 9a df 6b a6 fe 3f 98 66 09 3b 2f ed
0003c0 cc 69 18 0c 03 0e 69 9e 5e 84 fd df dd ab e6 f8
0003d0 6e 25 db 56 6a 6b 84 97 2b dd ee 7f b4 55 54 32
0003e0 00 18 f1 0a aa aa 30 a8 88 8a a8 a8 a0 05 55 55
0003f0 1e 00 f5 44 fb 2c e1 0a a3 ed b3 ff d9

Marker segment

JPEG Fileは0xFFではじまる2byteのMarkerによって区切られます。
Markerと"Markerで区切られるデータ"を併せて Marker segmentといいます。
次のMarkerがあります。

Name	Value
SOI (Start of image)	0xFFD8
APPn (Reserved for application segments)	0xFFE0 - 0xFFEF
DQT (Define quantization table)	0xFFDB
SOF0 (Start Of Frame markers;Baseline DCT)	0xFFC0
DHT (Define Huffman table)	0xFFC4
SOS (Start of scan)	0xFFDA
EOI (End of image)	0xFFD9

JPEG FileはSOIではじまり、EOIで終わります。Segmentの並びは次の通りです。
画像データ本体は ECS(entropy-coded segment) にあります。
ECSにはMarkerがありません。SOS直後に置かれます。

次にそれぞれのMarker segmentを説明します。

SOF0 (Start Of Frame marker;Baseline DCT)

このSegmentには次の情報があります。

Parameter	Size (bits)	Values
Marker	16	0xFFC0
Lf	16	Lfも含むlength. 8+3*Nf
P	8	the precision in bits for the samples
Y	16	Number of lines
X	16	Number of samples per line
Nf	8	Number of image components
Ci	8	Component identifier
Hi	4	Horizontal sampling factor
Vi	4	Vertical sampling factor
Tqi	8	Quantization table destination selector

earth.jpg の値は次の通りです。

ff c0 00 11 08 00 1f 00 1e 03 01 22 00 02 11 01 03 11 01

SOF:
  Lf=17
  P =8
  Y =31
  X =30
  Nf=3
    C(1) H 2 V 2 Tq 0 
    C(2) H 1 V 1 Tq 1 
    C(3) H 1 V 1 Tq 1 
  Hmax=2 Vmax=2

SOS (Start of scan)

このSegmentには次の情報があります。
SOS は各 component が使用する entropy coding tableを指定します。

Parameter	Size (bits)	Values
Marker	16	0xFFDA
Ls	16	Lsも含むlength. 6+2*Ns
Ns	8	Number of image components
Csj	8	Scan component selector
Tdj	4	DC entropy coding table destination selector
Taj	4	AC entropy coding table destination selector
Ss	8	Start of spectral or predictor selection
Se	8	End of spectral selection
Ah	4	Successive approximation bit position high
Al	4	Successive approximation bit position low or point transform

earth.jpg の値は次の通りです。

ff da 00 0c 03 01 00 02 11 03 11 00 3f 00

SOS:
  Ls=12
  Ns=3
  Ss=0
  Se=63
  C(1) Td=0 Ta=0
  C(2) Td=1 Ta=1
  C(3) Td=1 Ta=1
  Ah=0
  Al=0

DQT (Define quantization table)

このSegmentには次の情報があります。

Parameter	Size (bits)	Values
Marker	16	0xFFDB
Lq	16	Lqも含むlength. 67
Pq	4	Quantization table element precision.0:8bit,1:16bit
Tq	4	Quantization table destination identifier
Qk	8	Quantization table element

earth.jpg の値は次の通りです。2つのDQTがあります。
それぞれ Y'、Cb/Cr の Quantization tableです。

ff db 00 43 00 02 01 01 02 01 01 02 02 02 02 02 02 02 02 03 05 03 03 03 03 03 06 04 04 03 05 07 06 07 07 07 06 07 07 08 09 0b 09 08 08 0a 08 07 07 0a 0d 0a 0a 0b 0c 0c 0c 0c 07 09 0e 0f 0d 0c 0e 0b 0c 0c 0c
ff db 00 43 01 02 02 02 03 03 03 06 03 03 06 0c 08 07 08 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c 0c

DQT:
  Lq=67
  Pq=0(8bit)
  Tq=0
  Q
  002 001 001 002 001 001 002 002 002 002 002 002 002 002 003 005 003 003 003 003 003 006 004 004 003 005 007 006 007 007 007 006 007 007 008 009 011 009 008 008 010 008 007 007 010 013 010 010 011 012 012 012 012 007 009 014 015 013 012 014 011 012 012 012 

DQT:
  Lq=67
  Pq=0(8bit)
  Tq=1
  Q
  002 002 002 003 003 003 006 003 003 006 012 008 007 008 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012

DHT (Define Huffman table)

このSegmentには次の情報があります。

Parameter	Size (bits)	Values
Marker	16	0xFFC4
Lh	16	Lhも含むlength.
Tc	4	0 = DC table, 1 = AC table
Th	4	Huffman table destination identifier
Li	8	Number of Huffman codes of length i
Vi,j	8	Value associated with each Huffman code

earth.jpg の値は次の通りです。

ff c4 00 17 00 01 01 01 01 00 00 00 00 00 00 00 00 00 00 00 00 09 08 0a 05
ff c4 00 2f 10 00 02 01 03 03 03 03 01 07 05 00 00 00 00 00 00 01 02 03 04 05 11 06 07 12 08 09 21 00 13 41 0a 14 15 22 23 31 32 42 51 52 61 71 91
ff c4 00 17 01 00 03 01 00 00 00 00 00 00 00 00 00 00 00 00 00 00 02 04 03
ff c4 00 1a 11 00 02 03 01 01 00 00 00 00 00 00 00 00 00 00 00 00 01 02 11 12 03 21

DHT:
  Lh=23
  Tc=0(DC)
  Th=0
  L(BITS)
  001 001 001 001 000 000 000 000 000 000 000 000 000 000 000 000 
  HUFFSIZE last=4
  001 002 003 004 
  HUFFVAL(V) HUFFCODE
  09         0
  08         10
  0A         110
  05         1110
  MINCODE
  0 2 6 14 0 0 0 0 0 0 0 0 0 0 0 0 
  MAXCODE
  0 2 6 14 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 
  VALPTR
  0 1 2 3 0 0 0 0 0 0 0 0 0 0 0 0 

DHT:
  Lh=47
  Tc=1(AC)
  Th=0
  L(BITS)
  000 002 001 003 003 003 003 001 007 005 000 000 000 000 000 000 
  HUFFSIZE last=28
  002 002 003 004 004 004 005 005 005 006 006 006 007 007 007 008 
  009 009 009 009 009 009 009 010 010 010 010 010 
  HUFFVAL(V) HUFFCODE
  01         00
  02         01
  03         100
  04         1010
  05         1011
  11         1100
  06         11010
  07         11011
  12         11100
  08         111010
  09         111011
  21         111100
  00         1111010
  13         1111011
  41         1111100
  0A         11111010
  14         111110110
  15         111110111
  22         111111000
  23         111111001
  31         111111010
  32         111111011
  42         111111100
  51         1111111010
  52         1111111011
  61         1111111100
  71         1111111101
  91         1111111110
  MINCODE
  0 0 4 10 26 58 122 250 502 1018 0 0 0 0 0 0 
  MAXCODE
  -1 1 4 12 28 60 124 250 508 1022 -1 -1 -1 -1 -1 -1 
  VALPTR
  0 0 2 3 6 9 12 15 16 23 0 0 0 0 0 0 

DHT:
  Lh=23
  Tc=0(DC)
  Th=1
  L(BITS)
  000 003 001 000 000 000 000 000 000 000 000 000 000 000 000 000 
  HUFFSIZE last=4
  002 002 002 003 
  HUFFVAL(V) HUFFCODE
  00         00
  02         01
  04         10
  03         110
  MINCODE
  0 0 6 0 0 0 0 0 0 0 0 0 0 0 0 0 
  MAXCODE
  -1 2 6 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 
  VALPTR
  0 0 3 0 0 0 0 0 0 0 0 0 0 0 0 0 

DHT:
  Lh=26
  Tc=1(AC)
  Th=1
  L(BITS)
  000 002 003 001 001 000 000 000 000 000 000 000 000 000 000 000 
  HUFFSIZE last=7
  002 002 003 003 003 004 005 
  HUFFVAL(V) HUFFCODE
  00         00
  01         01
  02         100
  11         101
  12         110
  03         1110
  21         11110
  MINCODE
  0 0 4 14 30 0 0 0 0 0 0 0 0 0 0 0 
  MAXCODE
  -1 1 6 14 30 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 
  VALPTR
  0 0 2 5 6 0 0 0 0 0 0 0 0 0 0 0

Huffman table 作成

Li と Vi,j から Huffman tableを作成する方法を説明します。

Li は length i の Huffman code の個数を表します。iの範囲は[1, 16]です。
Vi,j は length i の Huffman code の j 番目の値です。

Li の具体値が次のように与えられているとします。

L(BITS)
000 002 003 001 001 000 000 000 000 000 000 000 000 000 000 000

1bit の code が 0個
2bit の code が 2個
3bit の code が 3個
4bit の code が 1個
5bit の code が 1個
あることが分かります。合計7個です。

次に HUFFSIZE を作成します。HUFFSIZE はcodeのbit数を左詰めで配列にしたものです。
個数は7個です。

HUFFSIZE last=7
002 002 003 003 003 004 005

次の手順で HUFFSIZE から Huffman code を計算します。

code = 0

HUFFSIZE[0]の値を取り出します。 2(bit)です。
code を LSB から 2bit 取り出して Huffman code 00b とします。
code++します。

HUFFSIZE[1]の値を取り出します。 2(bit)です。
code を LSB から 2bit 取り出して Huffman code 01b とします。
code++します。

HUFFSIZE[2]の値を取り出します。 3(bit)です。
bit数が2から3に増えたため codeを 1bit 左シフト(code<<1)します。
code を LSB から 3bit 取り出して Huffman code 100b とします。
code++します。

上記をHUFFSIZE[6]まで繰り返します。

得られたHuffman codeをHUFFCODEと表記します。
HUFFCODEに対応する値が Vi,j です。 HUFFVALと表記します。

  HUFFVAL(V) HUFFCODE
  00         00
  01         01
  02         100
  11         101
  12         110
  03         1110
  21         11110

ECS

ECS は Entropy coding した MCU の並びです。
MCU は data unit の並びです。
data unit は zig-zag order に並んでいます。ZZ(0) ~ ZZ(63)と表記します。
ZZ(0)をDC, ZZ(1) ~ ZZ(63) AC といいます。
DC と AC は異なる Entropy codingを行います。次にそれぞれのEntropy codingを説明します。

DC

DCは直前のdata unitのDCからの差分を符号化します。直前のDCをPREDとします。
次の式でZZ(0)を求めます。DIFFが符号化されるデータそのものです。PREDの初期値は0です。
PREDはcomponent毎にあります。

ZZ(0) = PRED + DIFF

EncoderはDIFFの範囲に応じてSSSSを決めます。対応は次の通りです。

SSSS	DIFF
0	0
1	-1,1
2	-3,-2,2,3
3	-7..-4,4..7
4	-15..-8,8..15
5	-31..-16,16..31
6	-63..-32,32..63
7	-127..-64,64..127
8	-255..-128,128..255
9	-511..-256,256..511
10	-1 023..-512,512..1 023
11	-2 047..-1 024,1 024..2 047

SSSSをHuffman codingします。
SSSSのHuffman code の後ろに値を表すbit sequence(DIFF'とします) を SSSS bit 並べます。
Huffman code と SSSSの対応表 (Huffman table) は DHT segmentに配置されます。

具体例で説明します。
Huffman table と bit sequenceが次のように与えられているとします。
ZZ(0)を求めます。

Huffman table
  HUFFVAL(V) HUFFCODE
  09         0
  08         10
  0A         110
  05         1110

bit sequence
  10 00001101

bit sequence の先頭は 10b です。
したがってHuffman tableからSSSSは8であることが分かります。
次の 8bit を取り出します。00001101b です。10進数で13です。
SSSS が 8 の場合、DIFFの範囲は -255..-128,128..255 です。
DIFF'がそのままDIFFとなるわけではありません。
次の式により、2^(SSSS-1) = 2^7 = 128より小さい場合は、負の数を表します。

DIFF' = (Huffman code の後ろにある 値を表すbit sequence);
VT = 2^(SSSS-1);
if (DIFF' < VT) {
  DIFF = DIFF' + (-1 << T) + 1;
} else {
  DIFF = DIFF';
}

DIFF = 13 + (-2^8) + 1 = 13 - 256 + 1 = -242
DIFFは-242になります。DIFFとPRED から ZZ(0)を求めます。その後、PRED = ZZ(0)とします。

AC

AC では 0が続く頻度が高いため、

"0の数"(上位4bit) RRRRとします。
"0の後ろにある値のbit length"(下位4bit) SSSSとします。

を　Huffman coding します。SSSSについてはDCと同じ処理を行います。
Huffman code と RRRR SSSSの対応表 (Huffman table) は DHT segmentに配置されます。

EOB(end-of-block)は以降のZZが0であることを示します。data unitの終端です。値は0x00です。
ZRLは0が16個あることを表します。値は0xF0です。

具体例で説明します。
Huffman table と bit sequenceが次のように与えられているとします。
ZZ(1) ~ ZZ(63)を求めます。

Huffman table
  HUFFVAL(V) HUFFCODE
  01         00
  02         01
  03         100
  04         1010
  05         1011
  11         1100
  06         11010
  07         11011
  12         11100
  08         111010
  09         111011
  21         111100
  00         1111010
  13         1111011
  41         1111100
  0A         11111010
  14         111110110
  15         111110111
  22         111111000
  23         111111001
  31         111111010
  32         111111011
  42         111111100
  51         1111111010
  52         1111111011
  61         1111111100
  71         1111111101
  91         1111111110

bit sequence
  111111001 011 1111010

bit sequence の先頭は 111111001b です。
Huffman tableからRRRR=2 SSSS=3であることが分かります。
次の 3bit を取り出します。011b です。10進数で3です。
SSSS が 3 の場合、DIFFの範囲は -7..-4,4..7 です。
2^(SSSS-1) = 2^2 = 4 より小さい場合は、負の数を表します。
DIFF = 3 + (-2^3) + 1 = 3 - 8 + 1 = -4
DIFFは-4になります。RRRR=2により0が2つ続きます。ZZ(1) = ZZ(2) = 0です。
ZZ(3) = DIFF = -4になります。
次のbit sequenceは1111010です。
Huffman tableからEOBあることが分かります。
よってZZ(4)~ZZ(63)は0になります。

0xFF00の扱い

Entropy codeとMarkerが同じ値になる可能性があります。
これを避けるためにEncoderは0xFFが現れる Entropy code の後ろに0x00を挿入します。
Decoderは0xFF00を検出すると後ろの0x00を除きます。

Decoding process

次の順番でdecodeします。

Marker detection
Entropy decoding
Dequantization
IDCT
Level shift

Marker detection

FileからMarker segmentを取り出します。

SOF0 : X, Y, Hi, Vi, Tqi を取得します。i=1,2,3です。
SOS : Tdi, Tai を取得します。i=1,2,3です。
DQT : Y' / CbCr それぞれの Quantization table を取得します。
DHT : Y' DC / Y' AC / CbCr DC / CbCr AC それぞれのHuffman tableを取得します。

SOS 直後からECSがはじまります。

Entropy decoding

ECSをEntropy decodeします。

Dequantization

Entropy decode した data に対して Quantization tableを掛けます。

IDCT

Dequantize 後のdata を zig-zag order に沿って 8x8 blockに戻します。
8x8 block 対してIDCTを行います。
計算方法は以下を参照ください。
https://qiita.com/tobira-code/items/91f3578cd7ed5b19c1f9

Level shift

Encoderはimage data を -128 level shiftしたものをFDCTします。
DecoderはIDCT後のdata を +128 level shiftします。

sample code

sample codeを示します。実行環境はMINGW64です。

$ gcc --version
gcc.exe (Rev2, Built by MSYS2 project) 6.2.0

makefile

CFLAGS=-I. -Wall -g -mno-ms-bitfields
INCS=
OBJS=test.o
LIBS=
TARGET=test

all: $(TARGET)

%.o: %.c $(INCS)
	$(CC) $(CFLAGS) -c -o $@ $<

$(TARGET): $(OBJS)
	$(CC) $(CFLAGS) -o $@ $^ $(LIBS)

clean:
	rm -rf $(TARGET) *.o

test.c

# include <math.h>
# include <stdbool.h>
# include <stdint.h>
# include <stdio.h>
# include <stdlib.h>
# include <string.h>

/* Marker code */
# define SOF0 (0xC0) /* Baseline DCT                         */
# define DHT  (0xC4) /* Define Huffman table(s)              */
# define SOS  (0xDA) /* Start of scan                        */
# define DQT  (0xDB) /* Define quantization table(s)         */

/* Frame header */
struct sofn_c_t {
  uint8_t C;    /* Component identifier                     */
  uint8_t V:4;  /* Vertical sampling factor                 */
  uint8_t H:4;  /* Horizontal sampling factor               */
  uint8_t Tq;   /* Quantization table destination selector  */
} __attribute__((packed));

struct sofn_t {
  uint16_t Lf;  /* Frame header length                      */
  uint8_t P;    /* Sample precision                         */
  uint16_t Y;   /* Number of lines                          */
  uint16_t X;   /* Number of samples per line               */
  uint8_t Nf;   /* Number of image components in frame      */
# define MAX_NF (3)
  struct sofn_c_t C[MAX_NF];
  uint8_t Hmax;
  uint8_t Vmax;
} __attribute__((packed));

/* Scan header */
struct sos_c_t {
  uint8_t Cs;   /* Scan component selector                      */
  uint8_t Ta:4; /* AC entropy coding table destination selector */
  uint8_t Td:4; /* DC entropy coding table destination selector */
} __attribute__((packed));

struct sos_t {
  uint16_t Ls;  /* Scan header length                           */
  uint8_t Ns;   /* Number of image components in scan           */
# define MAX_NS (3)
  struct sos_c_t C[MAX_NS];
  uint8_t Ss;   /* Start of spectral or predictor selection     */
  uint8_t Se;   /* End of spectral selection                    */
  uint8_t Al:4; /* Successive approximation bit position low or point transform */
  uint8_t Ah:4; /* Successive approximation bit position high                   */
} __attribute__((packed));

/* Quantization table */
# define PQ_8BIT  (0)
struct dqt_t {
  uint16_t Lq;  /* Quantization table definition length       */
  uint8_t Tq:4; /* Quantization table destination identifier  */
  uint8_t Pq:4; /* Quantization table element precision       */
# define Q_LEN (128)
  uint8_t Q[Q_LEN];
} __attribute__((packed));

/* Huffman table */
struct dht_t {
  uint16_t Lh;      /* Huffman table definition length                              */
  uint8_t Th:4;     /* Huffman table destination identifier                         */
  uint8_t Tc:4;     /* Table class - 0 = DC table or lossless table, 1 = AC table.  */
# define L_LEN (16)
  uint8_t L[L_LEN]; /* Number of Huffman codes of length i                          */
# define V_LEN (256)
  uint8_t V[V_LEN]; /* Value associated with each Huffman code                      */
  uint8_t HUFFSIZE[V_LEN];
  uint8_t last;
  uint16_t HUFFCODE[V_LEN];
  uint8_t VALPTR[V_LEN];
  int MINCODE[V_LEN];
  int MAXCODE[V_LEN];
} __attribute__((packed));

struct info_c_t {
  int dc_pred;
  uint8_t* component;
  int width;
  int height;
} __attribute__((packed));

struct info_t {
  struct info_c_t C[MAX_NF];
} __attribute__((packed));

static FILE* g_fp;
static struct info_t g_info;

static int decode(struct dht_t* dht, int* code, int* size);
static int extend(int V, int T);
static int getbytes(uint8_t* buf, int size);
static int nextbit();
static int receive(int ssss);

static const int zig_zag_order[][8] = {
  {  0,  1,  5,  6, 14, 15, 27, 28 },
  {  2,  4,  7, 13, 16, 26, 29, 42 },
  {  3,  8, 12, 17, 25, 30, 41, 43 },
  {  9, 11, 18, 24, 31, 40, 44, 53 },
  { 10, 19, 23, 32, 39, 45, 52, 54 },
  { 20, 22, 33, 38, 46, 51, 55, 60 },
  { 21, 34, 37, 47, 50, 56, 59, 61 },
  { 35, 36, 48, 49, 57, 58, 62, 63 },
};

static void IDCT(int F[][8], int f[][8])
{
  int x,y,u,v;
  double N=8.0;
  
  for (y=0; y<8; y++) {
    for (x=0; x<8; x++) {
      double val = 0;
      for (v=0; v<8; v++) {
        for (u=0; u<8; u++) {
          double cu,cv;
          if (u == 0) {
            cu = 0.707106781; // 1/2^(1/2)
          } else {
            cu = 1.0;
          }
          if (v == 0) {
            cv = 0.707106781; // 1/2^(1/2)
          } else {
            cv = 1.0;
          }
          val += cu*cv*F[v][u]*cos((2.0*x+1.0)*(double)u*M_PI/(2.0*N))*cos((2.0*y+1.0)*(double)v*M_PI/(2.0*N));
        }
      }
      f[y][x] = round((2.0/N)*val);
    }
  }
}

static void print_binary(uint16_t value, int len)
{
  int i;
  
  for (i=0; i<len && i<16; i++) {
    int shift = len - i - 1;
    printf("%d", (value>>shift)&1);
  }
}

static void swap2byte(uint8_t* pbuf)
{
  uint8_t tmp;
  tmp = pbuf[1];
  pbuf[1] = pbuf[0];
  pbuf[0] = tmp;
}

static int getbytes(uint8_t* buf, int size)
{
  size_t sz = fread(buf, 1, size, g_fp);
  if (sz != size) {
    if (feof(g_fp)) {
      return 1; /* EOF */
    } else {
      return -1;
    }
  }
  
  return 0;
}

static void print_frame_header(struct sofn_t* sofn)
{
  int i;
  
  printf("SOF:\n");
  printf("  Lf=%d\n", sofn->Lf);
  printf("  P =%d\n", sofn->P);
  printf("  Y =%u\n", sofn->Y);
  printf("  X =%u\n", sofn->X);
  printf("  Nf=%u\n", sofn->Nf);
  for (i=0; i<sofn->Nf; i++) {
    printf("    C(%d) H %d V %d Tq %d \n", sofn->C[i].C, sofn->C[i].H, sofn->C[i].V, sofn->C[i].Tq);
  }
  printf("  Hmax=%u Vmax=%u\n", sofn->Hmax, sofn->Vmax);
  printf("\n");
}

static bool is_known_marker(uint8_t marker)
{
  if ((marker == SOF0) || (marker == DHT ) || (marker == SOS ) || (marker == DQT )) {
    return true;
  } else {
    return false;
  }
}

static uint16_t next_marker()
{
  int ret;
  uint16_t marker;
  
  while (1) {
    ret = getbytes((uint8_t*)&marker, 2);
    if (ret != 0) {
      return ret;
    }
    if ((marker & 0x00FF) == 0xFF && (is_known_marker(marker>>8))) {
      break;
    }
  }
  
  return marker>>8;
}

static void parse_frame_header(struct sofn_t* sofn)
{
  int ret, offset = 0, i;
  uint16_t marker;
  uint8_t *pbuf = (uint8_t*)sofn;
  
  marker = next_marker();
  if (marker != SOF0) {
    printf("error parse_frame_header %d\n", marker);
    return;
  }
  
  memset(pbuf, 0, sizeof(struct sofn_t));
  
  ret = getbytes(pbuf+offset, 2);
  if (ret != 0) {
    return;
  }
  swap2byte(pbuf+offset);
  offset += 2;
  
  ret = getbytes(pbuf+offset, sofn->Lf - 2);
  if (ret != 0) {
    return;
  }
  
  /* P */
  offset += 1;
  /* Y */
  swap2byte(pbuf+offset);
  offset += 2;
  /* X */
  swap2byte(pbuf+offset);
  offset += 2;
  /* Nf */
  offset += 1;
  /* calc Hmax Vmax */
  for (i=0; i<sofn->Nf && i<MAX_NF; i++) {
    if (sofn->Hmax < sofn->C[i].H) {
      sofn->Hmax = sofn->C[i].H;
    }
    if (sofn->Vmax < sofn->C[i].V) {
      sofn->Vmax = sofn->C[i].V;
    }
  }
}

static void print_quantization_table(struct dqt_t* dqt)
{
  int i;
  
  printf("DQT:\n");
  printf("  Lq=%d\n", dqt->Lq);
  printf("  Pq=%d(%s)\n", dqt->Pq, dqt->Pq == 0 ? "8bit" : "16bit");
  printf("  Tq=%d\n", dqt->Tq);
  printf("  Q\n  ");
  if (dqt->Pq == PQ_8BIT) {
    for (i=0; i<Q_LEN/2; i++) {
      printf("%03d ", dqt->Q[i]);
    }
  }
  printf("\n\n");
}

static void parse_quantization_table(struct dqt_t* dqt)
{
  int ret, offset = 0;
  uint16_t marker;
  uint8_t *pbuf = (uint8_t*)dqt;
  
  marker = next_marker();
  if (marker != DQT) {
    printf("error parse_quantization_table %d\n", marker);
  }
  
  memset(pbuf, 0, sizeof(struct dqt_t));
  
  ret = getbytes(pbuf+offset, 2);
  if (ret != 0) {
    return;
  }
  swap2byte(pbuf+offset);
  offset += 2;
  
  ret = getbytes(pbuf+offset, dqt->Lq - 2);
  if (ret != 0) {
    return;
  }
  offset += dqt->Lq - 2;
}

static void print_scan_header(struct sos_t* sos)
{
  int i;
  
  printf("SOS:\n");
  printf("  Ls=%d\n", sos->Ls);
  printf("  Ns=%d\n", sos->Ns);
  printf("  Ss=%d\n", sos->Ss);
  printf("  Se=%d\n", sos->Se);
  for (i=0; i<sos->Ns && i<MAX_NS; i++) {
    printf("  C(%d) Td=%d Ta=%d\n", sos->C[i].Cs, sos->C[i].Td, sos->C[i].Ta);
  }
  printf("  Ah=%d\n", sos->Ah);
  printf("  Al=%d\n", sos->Al);
  printf("\n");
}

static void generate_size_table(struct dht_t* dht)
{
  int k = 0 , i = 1, j = 1;
  
  while (1) {
    if (dht->L[i-1] < j && i-1 < L_LEN) {
      i = i + 1;
      j = 1;
      
      if (16 < i) {
        dht->HUFFSIZE[k] = 0;
        dht->last = k;
        break;
      }
    } else if (k < V_LEN) {
      dht->HUFFSIZE[k] = i;
      k = k + 1;
      j = j + 1;
    } else {
      printf("error fatal generate_size_table\n");
      break; /* error */
    }
  }
}

static void generate_code_table(struct dht_t* dht)
{
  int k = 0, code = 0, si = dht->HUFFSIZE[0];
  
  while (1) {
    dht->HUFFCODE[k] = code;
    code = code + 1;
    k = k + 1;
    
    if (dht->HUFFSIZE[k] != si) {
      if (dht->HUFFSIZE[k] == 0) {
        break;
      } else {
        do {
          code = code << 1;
          si = si + 1;
        } while (dht->HUFFSIZE[k] != si);
      }
    }
  }
}

static int decoder_tables(struct dht_t* dht)
{
  int i = 0, j = 0;
  
  for (i=0; i<16; i++) {
    if (dht->L[i] == 0) {
      dht->MAXCODE[i] = -1;
    } else {
      dht->VALPTR[i] = j;
      dht->MINCODE[i] = dht->HUFFCODE[j];
      j = j + dht->L[i] - 1;
      dht->MAXCODE[i] = dht->HUFFCODE[j];
      j = j + 1;
    }
  }
  
  return 0;
}

static void parse_huffman_table(struct dht_t* dht)
{
  int ret, offset = 0;
  uint16_t marker;
  uint8_t *pbuf = (uint8_t*)dht;
  
  marker = next_marker();
  if (marker != DHT) {
    printf("error parse_huffman_table %d\n", marker);
    return;
  }
  
  memset(pbuf, 0, sizeof(struct dht_t));
  
  ret = getbytes(pbuf+offset, 2);
  if (ret != 0) {
    return;
  }
  swap2byte(pbuf+offset);
  offset += 2;
  
  ret = getbytes(pbuf+offset, dht->Lh - 2);
  if (ret != 0) {
    return;
  }
  offset += dht->Lh - 2;
  
  generate_size_table(dht);
  generate_code_table(dht);
  decoder_tables(dht);
}

static void print_huffman_table(struct dht_t* dht)
{
  int i, sum = 0;
  
  printf("DHT:\n");
  printf("  Lh=%d\n", dht->Lh);
  printf("  Tc=%d(%s)\n", dht->Tc, dht->Tc==0 ? "DC" : "AC");
  printf("  Th=%d\n", dht->Th);
  printf("  L(BITS)\n  ");
  for (i=0; i<L_LEN; i++) {
    sum += dht->L[i];
    printf("%03d ", dht->L[i]);
  }
  printf("\n");
  
  printf("  HUFFSIZE last=%d\n  ", dht->last);
  for (i=0; i<dht->last; i++) {
    printf("%03d ", dht->HUFFSIZE[i]);
    if ((i+1)%16 == 0 && (i+1) != dht->last)
      printf("\n  ");
  }
  printf("\n");
  
  printf("  HUFFVAL(V) HUFFCODE\n");
  for (i=0; i<dht->last; i++) {
    printf("  %02X         ", dht->V[i]);
    print_binary(dht->HUFFCODE[i], dht->HUFFSIZE[i]);
    printf("\n");
  }
  
  printf("  MINCODE\n  ");
  for (i=0; i<16; i++) {
    printf("%d ", dht->MINCODE[i]);
  }
  printf("\n");
  
  printf("  MAXCODE\n  ");
  for (i=0; i<16; i++) {
    printf("%d ", dht->MAXCODE[i]);
  }
  printf("\n");
  
  printf("  VALPTR\n  ");
  for (i=0; i<16; i++) {
    printf("%d ", dht->VALPTR[i]);
  }
  printf("\n\n");
}

static void parse_scan_header(struct sos_t* sos)
{
  int ret, offset = 0;
  uint16_t marker;
  uint8_t *pbuf = (uint8_t*)sos;
  
  marker = next_marker();
  if (marker != SOS) {
    printf("error parse_scan_header %d\n", marker);
    return;
  }
  
  memset(pbuf, 0, sizeof(struct sos_t));
  
  ret = getbytes(pbuf+offset, 2);
  if (ret != 0) {
    return;
  }
  swap2byte(pbuf+offset);
  offset += 2;
  
  ret = getbytes(pbuf+offset, sos->Ls - 2);
  if (ret != 0) {
    return;
  }
  offset += sos->Ls - 2;
}

static void decode_ac_coefficients(struct dht_t* dht, int* ZZ)
{
  int K = 0, code, size;
  
  while (1) {
    int RS = decode(dht, &code, &size);
    int SSSS = RS % 16;
    int RRRR = RS >> 4;
    int R = RRRR;
    
    if (SSSS == 0) {
      if (R == 15) {
        K += 16;
        if (62 <= K) {
          break;
        }
      } else {
        break;
      }
    } else {
      K += R;
      ZZ[K] = receive(SSSS);
      ZZ[K] = extend(ZZ[K], SSSS);
      if (62 <= K) {
        break;
      } else {
        K += 1;
      }
    }
  }
}

static int decode_dc_coefficient(struct dht_t* dht)
{
  int code, size;
  int T = decode(dht, &code, &size);
  int DIFF = receive(T);
  return extend(DIFF, T);
}

static void decode_entropy_coded_data_unit(int index, struct dht_t* dc, struct dht_t* ac, int* ZZ)
{
  ZZ[0] = decode_dc_coefficient(dc) + g_info.C[index].dc_pred;
  g_info.C[index].dc_pred = ZZ[0];
  decode_ac_coefficients(ac, &ZZ[1]);
}

static uint8_t nextbyte()
{
  int ret;
  uint8_t val = 0;
  
  ret = getbytes(&val, 1);
  if (ret != 0) {
    exit(0);
  }
  return val;
}

static void level_shift(int f[8][8])
{
  int i,j;
  
  for (i=0; i<8; i++) {
    for (j=0; j<8; j++) {
      f[i][j] += 128;
      if (f[i][j] < 0) {
        f[i][j] = 0;
      }
      if (255 < f[i][j]) {
        f[i][j] = 255;
      }
    }
  }
}

static void rearrange_zig_zag_order(int* ZZ, int F[8][8])
{
  int i,j;
  
  for (i=0; i<8; i++) {
    for (j=0; j<8; j++) {
      F[i][j] = ZZ[zig_zag_order[i][j]];
    }
  }
}

static void dequantize(struct dqt_t* dqt, int* ZZ)
{
  int i;
  
  for (i=0; i<64; i++) {
    ZZ[i] *= dqt->Q[i];
  }
}

static void copy_data_unit(int dux, int duy, int index, int f[8][8])
{
  int i,j;
  struct info_c_t* C = g_info.C;
  int stride = C[index].width;
  uint8_t* component = C[index].component;
  int offset = duy * 8 * stride + dux * 8;
  
  for (i=0; i<8; i++) {
    for (j=0; j<8; j++) {
      uint8_t* p = component + offset + stride * i;
      p[j] = f[i][j];
    }
  }
}

static void decode_mcu_component(int mcux, int mcuy, int index, uint8_t H, uint8_t V, struct dht_t* dc, struct dht_t* ac, struct dqt_t* dqt)
{
  int row, col;
  int ZZ[64];
  int F[8][8];
  int f[8][8];
  
  for (row=0; row<V; row++) {
    for (col=0; col<H; col++) {
      memset(ZZ, 0, sizeof(ZZ));
      decode_entropy_coded_data_unit(index, dc, ac, ZZ);
      dequantize(dqt, ZZ);
      rearrange_zig_zag_order(ZZ, F);
      IDCT(F, f);
      level_shift(f);
      copy_data_unit(mcux*H+col, mcuy*V+row, index, f);
    }
  }
}

static struct dht_t* get_dht(int index, struct dht_t* dht, int dht_size, bool is_dc)
{
  int i;
  
  for (i=0; i<dht_size; i++) {
    if (is_dc) {
      if ((dht[i].Tc == 0 /* DC */) && (dht[i].Th == index)) {
        return &dht[i];
      }
    } else {
      if ((dht[i].Tc == 1 /* AC */) && (dht[i].Th == index)) {
        return &dht[i];
      }
    }
  }
  
  return NULL;
}

static struct dqt_t* get_dqt(int index, struct dqt_t* dqt, int dqt_size)
{
  int i;
  
  for (i=0; i<dqt_size; i++) {
    if (dqt[i].Tq == index) {
      return &dqt[i];
    }
  }
  
  return NULL;
}

static void decode_mcu(int mcux, int mcuy, struct sofn_t* sofn, struct sos_t* sos, struct dht_t* dht, int dht_size, struct dqt_t* dqt, int dqt_size)
{
  int i;
  struct dht_t* dht_dc;
  struct dht_t* dht_ac;
  struct dqt_t* _dqt;
  
  for (i=0; i<sofn->Nf; i++) {
    int Td = sos->C[i].Td;
    int Ta = sos->C[i].Ta;
    int Tq = sofn->C[i].Tq;
    dht_dc = get_dht(Td, dht, dht_size, true  /* dc */);
    dht_ac = get_dht(Ta, dht, dht_size, false /* ac */);
    _dqt = get_dqt(Tq, dqt, dqt_size);
    decode_mcu_component(mcux, mcuy, i, sofn->C[i].H, sofn->C[i].V, dht_dc, dht_ac, _dqt);
  }
}

static void print_component(int index)
{
  uint8_t* component = g_info.C[index].component;
  int w = g_info.C[index].width;
  int h = g_info.C[index].height;
  int i,j;
  
  printf("P2\n");
  printf("%d %d\n", w, h);
  printf("255\n");
  for (i=0; i<h; i++) {
    for (j=0; j<w; j++) {
      printf("%3d ", component[w*i+j]);
    }
    printf("\n");
  }
  printf("\n");
}

static void decode_entropy_coded_segment(struct sofn_t* sofn, struct sos_t* sos, struct dht_t* dht, int dht_size, struct dqt_t* dqt, int dqt_size)
{
  int row, col, i;
  int xf = sofn->Hmax<<3; /* x8 */
  int yf = sofn->Vmax<<3;
  int mcux = (sofn->X+xf-1)/xf;
  int mcuy = (sofn->Y+yf-1)/yf;
  int x = mcux*xf;
  int y = mcuy*yf;
  
  printf("X=%d Y=%d Hmax=%d Vmax=%d x=%d y=%d mcux=%d mcuy=%d xf=%d yf=%d\n\n",
    sofn->X, sofn->Y, sofn->Hmax, sofn->Vmax, x, y, mcux, mcuy, xf, yf);
  
  for (i=0; i<sofn->Nf; i++) {
    int w = mcux * sofn->C[i].H * 8;
    int h = mcuy * sofn->C[i].V * 8;
    g_info.C[i].component = malloc(w*h);
    g_info.C[i].width = w;
    g_info.C[i].height = h;
  }
  
  for (row=0; row < mcuy; row++) {
    for (col=0; col < mcux; col++) {
      decode_mcu(col, row, sofn, sos, dht, dht_size, dqt, dqt_size);
    }
  }
  
  for (i=0; i<sofn->Nf; i++) {
    print_component(i);
    free(g_info.C[i].component);
  }
}

static int decode(struct dht_t* dht, int* code, int* size)
{
  int i = 0, j = 0, value = 0;
  int _code = 0;
  
  for (i=0; i<16; i++) {
    _code = (_code << 1) + nextbit();
    if (_code <= dht->MAXCODE[i]) {
      break;
    }
  }
  
  j = dht->VALPTR[i];
  j = j + _code - dht->MINCODE[i];
  value = dht->V[j]; /* V = HUFFVAL */
  if (code != NULL)
    *code = dht->HUFFCODE[j];
  if (size != NULL)
    *size = dht->HUFFSIZE[j];
  
  return value;
}

static int nextbit()
{
  static int cnt = 0;
  static uint8_t b = 0;
  int bit = 0;
  
  if (cnt == 0) {
    b = nextbyte();
    cnt = 8;
    
    if (b == 0xFF) {
      nextbyte(); /* skip zero byte */
    }
  }
  
  bit = b >> 7;
  cnt = cnt - 1;
  b = b << 1;
  
  return bit;
}

static int receive(int ssss)
{
  int i = 0, v = 0;
  
  for (i=0; i<ssss; i++) {
    v = (v << 1) + nextbit();
  }
  
  return v;
}

static int extend(int V, int T)
{
  int Vt = (1 << (T-1));
  
  if (V < Vt) {
    Vt = (-1 << T) + 1;
    V = V + Vt;
  }
  
  return V;
}

int main(int argc, char* argv[])
{
  struct dqt_t dqt[2];
  struct sofn_t sofn;
  struct dht_t dht[4];
  struct sos_t sos;
  
  if (argc != 2) {
    printf("option error %d\n", argc);
    return -1;
  }
  
  g_fp = fopen(argv[1], "rb");
  if (g_fp == NULL) {
    printf("error fopen %s\n", argv[1]);
    return -1;
  }
  
  /* quantization_table */
  parse_quantization_table(&dqt[0]); print_quantization_table(&dqt[0]);
  parse_quantization_table(&dqt[1]); print_quantization_table(&dqt[1]);
  
  /* frame_header */
  parse_frame_header(&sofn); print_frame_header(&sofn);
  
  /* huffman_table */
  parse_huffman_table(&dht[0]); print_huffman_table(&dht[0]);
  parse_huffman_table(&dht[1]); print_huffman_table(&dht[1]);
  parse_huffman_table(&dht[2]); print_huffman_table(&dht[2]);
  parse_huffman_table(&dht[3]); print_huffman_table(&dht[3]);
  
  /* scan_header */
  parse_scan_header(&sos); print_scan_header(&sos);
  
  /* entropy-coded segment */
  decode_entropy_coded_segment(&sofn, &sos, dht, 4, dqt, 2);
  
  fclose(g_fp);
  
  return 0;
}

ビルドして実行します。

console

$ make && ./test.exe earth.jpg
cc -I. -Wall -g -mno-ms-bitfields -c -o test.o test.c
cc -I. -Wall -g -mno-ms-bitfields -o test test.o
DQT:
  Lq=67
  Pq=0(8bit)
  Tq=0
  Q
  002 001 001 002 001 001 002 002 002 002 002 002 002 002 003 005 003 003 003 003 003 006 004 004 003 005 007 006 007 007 007 006 007 007 008 009 011 009 008 008 010 008 007 007 010 013 010 010 011 012 012 012 012 007 009 014 015 013 012 014 011 012 012 012

DQT:
  Lq=67
  Pq=0(8bit)
  Tq=1
  Q
  002 002 002 003 003 003 006 003 003 006 012 008 007 008 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012 012

SOF:
  Lf=17
  P =8
  Y =31
  X =30
  Nf=3
    C(1) H 2 V 2 Tq 0
    C(2) H 1 V 1 Tq 1
    C(3) H 1 V 1 Tq 1
  Hmax=2 Vmax=2

DHT:
  Lh=23
  Tc=0(DC)
  Th=0
  L(BITS)
  001 001 001 001 000 000 000 000 000 000 000 000 000 000 000 000
  HUFFSIZE last=4
  001 002 003 004
  HUFFVAL(V) HUFFCODE
  09         0
  08         10
  0A         110
  05         1110
  MINCODE
  0 2 6 14 0 0 0 0 0 0 0 0 0 0 0 0
  MAXCODE
  0 2 6 14 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1
  VALPTR
  0 1 2 3 0 0 0 0 0 0 0 0 0 0 0 0

DHT:
  Lh=47
  Tc=1(AC)
  Th=0
  L(BITS)
  000 002 001 003 003 003 003 001 007 005 000 000 000 000 000 000
  HUFFSIZE last=28
  002 002 003 004 004 004 005 005 005 006 006 006 007 007 007 008
  009 009 009 009 009 009 009 010 010 010 010 010
  HUFFVAL(V) HUFFCODE
  01         00
  02         01
  03         100
  04         1010
  05         1011
  11         1100
  06         11010
  07         11011
  12         11100
  08         111010
  09         111011
  21         111100
  00         1111010
  13         1111011
  41         1111100
  0A         11111010
  14         111110110
  15         111110111
  22         111111000
  23         111111001
  31         111111010
  32         111111011
  42         111111100
  51         1111111010
  52         1111111011
  61         1111111100
  71         1111111101
  91         1111111110
  MINCODE
  0 0 4 10 26 58 122 250 502 1018 0 0 0 0 0 0
  MAXCODE
  -1 1 4 12 28 60 124 250 508 1022 -1 -1 -1 -1 -1 -1
  VALPTR
  0 0 2 3 6 9 12 15 16 23 0 0 0 0 0 0

DHT:
  Lh=23
  Tc=0(DC)
  Th=1
  L(BITS)
  000 003 001 000 000 000 000 000 000 000 000 000 000 000 000 000
  HUFFSIZE last=4
  002 002 002 003
  HUFFVAL(V) HUFFCODE
  00         00
  02         01
  04         10
  03         110
  MINCODE
  0 0 6 0 0 0 0 0 0 0 0 0 0 0 0 0
  MAXCODE
  -1 2 6 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1
  VALPTR
  0 0 3 0 0 0 0 0 0 0 0 0 0 0 0 0

DHT:
  Lh=26
  Tc=1(AC)
  Th=1
  L(BITS)
  000 002 003 001 001 000 000 000 000 000 000 000 000 000 000 000
  HUFFSIZE last=7
  002 002 003 003 003 004 005
  HUFFVAL(V) HUFFCODE
  00         00
  01         01
  02         100
  11         101
  12         110
  03         1110
  21         11110
  MINCODE
  0 0 4 14 30 0 0 0 0 0 0 0 0 0 0 0
  MAXCODE
  -1 1 6 14 30 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1
  VALPTR
  0 0 2 5 6 0 0 0 0 0 0 0 0 0 0 0

SOS:
  Ls=12
  Ns=3
  Ss=0
  Se=63
  C(1) Td=0 Ta=0
  C(2) Td=1 Ta=1
  C(3) Td=1 Ta=1
  Ah=0
  Al=0

X=30 Y=31 Hmax=2 Vmax=2 x=32 y=32 mcux=2 mcuy=2 xf=16 yf=16

P2
32 32
255
  2   3   7   4   3   3   3   0   0   1   1   1   1  22  51  50  22   0   2   0   3   0   2   7   0   0   4   0   2   0   0   2
  3   8   0   0   0   2   0   4   3  17 139 241 252 254 255 253 255 252 233 130   3   2   2   1   7   5   0   8   0   3   6   2
  9   1   6   6   4   0   2   4 181 253 255 244 247 239 241 241 239 240 238 250 255 171   0   1   1   1   1   7   4   5   1   5
  0   3   2   0   0   0  88 254 255 243 246 244 245 246 235 240 245 239 239 235 224 239 255  71   1   6   0   3   7   0   2   0
  5   1   3   7   1 166 253 245 246 248 246 249 246 243 244 238 241 239 228 230 231 219 213 245 117   4   8   8   2   2   4   3
  1   3   6   0 152 255 246 252 246 250 244 245 250 243 240 248 242 237 236 228 236 227 206 235 243 103   0   5   3   7   1   9
  4  12   0  73 255 239 244 242 246 245 251 248 244 248 246 240 244 239 234 229 230 227 214 221 239 233  39   0   6   0   4   0
  7   0   1 254 241 240 248 248 250 252 248 250 252 247 245 250 238 238 235 235 236 225 239 205 212 199 191   3   3   0   2   0
  2   2 150 254 240 249 246 251 248 254 253 251 250 250 251 245 241 244 238 236 234 229 250 228 213 201 207  70   3   4   0   5
  3   0 255 245 243 241 254 248 248 253 252 251 251 250 250 246 244 243 243 233 230 219 229 211 231 194 190 175   0   3   1   1
  0 105 255 237 242 246 240 253 249 251 250 252 252 249 249 245 245 237 236 232 225 231 214 209 204 229 189 202  20   0   7   0
  3 196 246 245 243 242 249 242 249 250 249 252 252 248 247 244 239 245 237 234 230 213 218 209 207 188 227 179  89   5   0   3
  0 254 237 238 239 244 246 245 248 249 247 250 250 245 244 241 238 242 236 237 221 221 202 207 191 181 225 188 129   0   3   1
  0 255 234 244 248 242 242 239 246 247 246 248 247 243 242 239 237 233 231 228 228 219 205 196 191 183 178 172 176   0   0   1
 22 254 234 244 237 241 234 252 243 246 245 245 243 240 240 237 236 233 233 228 215 216 199 200 182 180 169 169 182   3   0   2
 33 255 227 236 237 236 240 237 241 245 244 243 241 239 240 236 234 234 223 220 215 207 214 205 181 170 186 153 176   1   2   1
 15 252 231 237 234 238 239 236 234 241 244 243 240 234 239 232 233 229 220 217 208 199 202 208 191 166 155 159 172   0   0   0
  1 255 235 233 230 241 246 241 238 235 239 248 252 243 238 238 234 225 216 212 210 198 187 190 180 171 153 158 162   0   0   0
  0 246 241 236 230 233 233 239 231 233 233 233 237 249 253 254 247 236 211 207 211 204 177 165 197 149 156 161 123   0   3   0
  5 157 249 236 221 231 233 229 225 228 234 239 229 230 242 250 248 240 210 206 220 220 197 170 198 160 144 165  69   0   0   1
  0  58 251 243 225 232 231 247 230 221 222 240 242 242 255 250 233 233 185 214 200 231 178 190 190 182 146 163  18   0   4   0
  8   2 239 249 224 229 228 220 244 234 233 243 249 240 234 212 213 206 188 193 191 244 207 216 197 194 148 120   1   0   0   0
  0   1  75 255 240 223 230 222 240 235 239 232 244 237 232 233 221 206 181 178 184 181 141 185 205 136 168  34   4   4   0   3
  4   8   0 216 232 223 238 225 237 241 239 204 232 242 218 204 215 192 202 203 168 148 152 156 154 141 118   0   2   7   6   7
  3   8   1   3 247 210 230 234 190 229 232 230 246 212 218 239 205 218 199 183 172 184 201 168 142 164   0   9   4   0   3   1
  4   4   1   2  24 251 204 231 216 224 231 237 242 221 241 232 224 215 187 199 190 189 166 172 176  38   0   0   4   5   3   3
  1   4   4   6   9  30 237 231 232 231 230 221 222 232 221 215 221 219 208 191 199 186 156 169  43   0   8  12   2   0   1   3
  3   9   3  10   4   0   6 191 220 209 216 216 224 218 217 214 209 203 201 187 190 173 161  23   2  11   3   0   0   2   3   0
  2   7   5   4   9   6   6   1  62 210 221 208 209 200 199 197 197 185 192 183 179  98   4   0   7   0   5   5   6   6   2   8
  4   4  11   3   3  11   6   8   2   2  52 116 183 211 206 201 195 179 130  68   0   0   5   8   2  11   0   1   6   0   0   1
  6   1   0  10   7   0   2   6   0   2   2   3   0   3   4   0   5   0   0   0   0  11   0   9   0   4   4   0   4   0   5   0
  3   6   0   8   4   2   0   7   1   2   3   0   6   0   1   1   2   0   2   0   0   5   3   4   0   1   2   5   2   2   0   0

P2
16 16
255
128 128 128 128 128 128 128 128 128 128 128 128 128 128 128 128
128 128 128 128 128 128 128 128 128 128 128 128 128 128 128 128
128 128 128 128 128 128 128 128 128 128 129 130 130 129 129 128
128 128 128 128 128 128 128 128 128 129 130 131 131 131 130 129
128 128 128 128 128 128 128 128 128 129 131 133 133 132 130 129
128 128 128 128 128 128 128 128 127 129 132 135 135 134 131 130
128 128 128 128 128 128 128 128 127 130 133 136 136 135 132 130
128 128 128 128 128 128 128 128 127 130 134 136 137 135 132 130
128 128 128 129 129 128 128 128 130 131 134 135 135 133 131 129
128 128 129 129 129 129 128 128 130 132 134 135 135 133 131 129
128 129 129 129 129 129 129 129 131 132 134 135 134 132 130 129
128 129 129 130 130 130 130 129 132 133 133 133 133 131 130 129
128 128 129 130 130 130 130 130 132 132 132 132 131 130 129 129
127 128 129 129 130 130 130 130 131 131 130 130 129 129 129 129
127 127 128 129 129 129 129 129 130 130 129 128 128 128 128 129
126 127 128 128 129 129 129 129 130 129 128 127 127 128 128 129

P2
16 16
255
128 128 128 128 128 128 128 128 128 128 128 128 128 128 128 128
128 128 128 128 128 128 128 128 128 128 128 128 128 128 128 128
128 128 128 128 128 128 128 128 128 128 128 128 128 128 128 128
128 128 128 128 128 128 128 128 128 128 127 127 127 127 128 128
128 128 128 128 128 128 128 128 128 128 127 127 127 127 128 128
128 128 128 128 128 128 128 128 128 127 127 126 126 127 127 128
128 128 128 128 128 128 128 128 128 127 126 126 126 126 127 128
128 128 128 128 128 128 128 128 128 127 126 126 126 126 127 128
128 128 128 128 128 128 128 128 128 128 127 126 126 126 127 128
128 128 128 128 128 128 128 128 128 127 126 126 126 127 127 128
128 128 128 128 128 128 128 128 127 127 126 126 126 127 128 129
128 128 128 128 128 128 128 128 127 127 127 127 127 128 128 129
128 128 128 128 128 128 128 128 127 127 127 127 127 128 129 129
128 128 128 128 128 128 128 128 127 127 127 128 128 128 128 129
128 128 128 128 128 128 128 128 128 128 128 128 128 128 128 128
128 128 128 128 128 128 128 128 128 128 128 128 128 128 128 128