2012년 2월 23일 목요일

모스부호(Morse code)와 Tip



간만에 블로그를 올린다.

이번에는 직접 제작한 매트랩 프로그램과 알면 유용한 매트랩 명령어에 대하여 알아보자!



그러기에 앞서, 항상 느끼는 것이지만,,,

매트랩은 공대생들에게 쉽고 유용한 필수도구라고 생각이든다.

짧은코딩으로, 쉽게 목적을 달성할 수 있는지 미리 살필 수 있고, 원리를 파악하기도 쉽다.

하드웨어랑 연결하고 조작도 쉬운것 같은데,,, 아직은 잘 모르겠다.

분명한것은 공부해야할 것들이 무궁무진하게 많고,

쉽게 코딩해서 강력한 기능을 맘대로 사용할 수 있다는 점이다.

수학과 공학 지식의 축적으로 이루어진 거대한 탑...^^

너무 높고, 방대해서 한 사람의 노력으로는 다 이해하기 어려운,,,

아무튼, 멋진 공학의 기초 도구이가 될것이다!ㅎ



그럼 이제, K5 자동차 광고에 이용해서 고급스러움과 신비스러움을 더 해줬던

모스부호의 세계로 빠져봅시다~

아래 광고의 '뚜~ 뚜 뚜~ 뚜 뚜 뚜 뚜 뚜' 의 세련된 음은

모스부호로써  'K5' 를 연속해서 가르키고 있다.


 


아래코드는 자작한 모스부호생성, 모스부호분석하여 변환하는 프로그램이다.</><>

아직 wave 파일 분석은 일반적이지 못해 좀더 수정해 볼 생각이다.



table={'a' '.-';
'b' '-...';
'c' '-.-.';
'd' '-..';
'e' '.';
'f' '..-.';
'g' '--.';
'h' '....';
'i' '..';
'j' '.---';
'k' '-.-';
'l' '.-..';
'm' '--';
'n' '-.';
'o' '---';
'p' '.--.';
'q' '--.-';
'r' '.-.';
's' '...';
't' '-';
'u' '..-';
'v' '...-';
'w' '.--';
'x' '-..-';
'y' '-.--';
'z' '--..';

'1' '.----';
'2' '..---';
'3' '...--';
'4' '....-';
'5' '.....';
'6' '-....';
'7' '--...';
'8' '---..';
'9' '----.';
'0' '-----'};

code=[];
c=input('문자를 입력해주세요!\n', 's'); c=lower(c);
for i=1:length(c)
for j=1:length(table)
if strcmp(c(i), table{j,1}), code=[code table{j,2}]; end
end
code=[code ' '];
end

fs=8000;
l=0.1;
dot_t=0:1/fs:l;      %
dash_t=0:1/fs:3*l;   % -
dot=sin(2*pi*1812.5*dot_t);
dash=sin(2*pi*1812.5*dash_t);
space=zeros(1,length(dot));
 

code
morse=[];
for i=1:length(code)
if strcmp(code(i),'.')
morse=[morse dot];
elseif strcmp(code(i),'-')
morse=[morse dash];
else
morse=[morse space space];
end
morse=[morse space];
end
sound(morse,fs)
% wavwrite(morse,fs,16, c)


% morese=wavread(c)
ev=(abs(hilbert(morse))>0.5);
% s(t)=Re[s_hat(t)·exp(jwt)]=si(t)cos(wt)-sq(t)sin(wt)  ,  s_hat(t)=si(t)+j·sq(t)
evs=downsample(ev,fs*l/8);
re=[]; k=1;
while k~=length(evs)
if length(evs)>=k+23 & sum(evs(k:k+23))==24, re=[re '-']; k=k+24;
elseif length(evs)>=k+23 & sum(evs(k:k+23))==0, re=[re ' ']; k=k+24;
elseif length(evs)>=k+7 & sum(evs(k:k+7))==8, re=[re '.']; k=k+8;
end
k=k+1;
end

id=[1 strfind(re, ' ')+1];
rec=[]; i=1; j=1;
for i=1:length(id)-1
for j=1:length(table)
if strcmp(re(id(i):id(i+1)-2), table{j,2}), rec=[rec table{j,1}]; end
end
end
upper(rec)



별건아니지만, hello world의 순수함으로 만들고, 뿌듯함을 느낀다.ㅎㅎ

다음은 잘 몰랐지만, 알고나면 굉장히 유용한 matlab의 주요한 기능을 살펴보자!



%-----------------------------------------------
% Create a new AVI file
%-----------------------------------------------
   t = linspace(0,2.5*pi,40);
       fact = 10*sin(t);
       fig=figure;
       aviobj = avifile('example.avi')
       [x,y,z] = peaks;
       for k=1:length(fact)
           h = surf(x,y,fact(k)*z);
           axis([-3 3 -3 3 -80 80])
           axis off
           caxis([-90 90])
           F = getframe(fig);
           aviobj = addframe(aviobj,F);
       end
       close(fig)
       aviobj = close(aviobj);
%-----------------------------------------------

%-----------------------------------------------
% Display wait bar
%-----------------------------------------------
h = waitbar(0,'Please wait...');
        for i=1:1000,
            % computation here %
            waitbar(i/1000,h)
        end
%-----------------------------------------------


 
%-----------------------------------------------
% Display legend
%-----------------------------------------------
x = 0:.2:12;
        plot(x,bessel(1,x),x,bessel(2,x),x,bessel(3,x));
        legend('First','Second','Third');
        legend('First','Second','Third','Location','NorthEastOutside')
%-----------------------------------------------


 

%-----------------------------------------------
% meshgrid, surf, contour
%-----------------------------------------------

[X,Y] = meshgrid(-1:.1:1);
        Z = 1-(X.^2 + Y.^2);
        surf(X,Y,Z)
        surfc(X,Y,Z)
[C, h]=contour(X, Y, Z);
clabel(C, h)
save('f.mat', 'X', 'Y', 'Z');
load('f.mat')
%-----------------------------------------------

 
%-----------------------------------------------
% Record and play back 2 seconds of single audio sampled at 8 kHz
%-----------------------------------------------
        Fs = 8000;
        y  = wavrecord(2*Fs, Fs, 'Single');
        wavplay(y, Fs);
%-----------------------------------------------


%-----------------------------------------------
% Calendar
%-----------------------------------------------
calendar(date)
%-----------------------------------------------

 
%-----------------------------------------------
% Capture
%-----------------------------------------------
vid = videoinput('winvideo', 1, 'YUY2_160x120');
src = getselectedsource(vid);
vid.FramesPerTrigger = 1;
start(vid);
vid.ReturnedColorspace = 'rgb';
nicelcw = getdata(vid);
imshow(nicelcw);
%-----------------------------------------------

 
 
%-----------------------------------------------
% Font color
%-----------------------------------------------
a='test';
text(1,1,a,'color','r')
%-----------------------------------------------




이상,,, 블로그에 업로드하기까지 오랜기간이 걸렸다. 그동안 항상 맘에 있었지만,,,

게을렀다고 할까?!ㅎ

한편으로는 허접한 블로그를 업로드하기위해, 나는 지난날을 그렇게 끙끙댔나보다~



올해,,, 만학도로 대학원을 입학했다.

몇일해본 대학원은 해주는 것을 받아먹는게 아니라,

노력하는 만큼 얻어가는 곳이라는 느낌이 확 든다.

공부해서,,, 졸업 쯔음엔 떳떳하고 멋있게 내 논문 발표할 수 있으면 좋겠다^^

그리고 좋은 블로그를 게시하여,

디지털통신 전공생의 가려운 부분을 박박 긁어주고 싶고, 지식도 나누고 싶다.

나 자신을 위하여 파이팅~!!ㅎ

2011년 11월 9일 수요일

채널코딩

채널코딩은 잡음에 오염된 채널을 통해 신뢰할 만한 통신을 성취하기 위하여, 소스 데이터를 변환하는 것을 말한다. 신뢰도를 높이는 대신 전송률과 대역폭, 복잡도는 절충해야 한다.

전송 오류제어하는 방법은,,,
수신측에서 독자적으로 오류를 일으킨 비트를 정정(FEC, 순방향 오류 정정),
수신단에서 오류 검출을 하고, 오류 검출될 때마다 역방향 링크를 통해 송신측에 보고하여 데이터를 재전송하도록 하는 방식(ARQ, 자동 반복 요구)

FEC 채널코딩을 분류하자면,,,
wave coding , block coding, convolution coding 이 있다.
블록 코딩은 현재의 k비트 메시지어는 과거의 k비트 메시지어와 관계 없이 독립적으로 출력을 만드어 내는데 반하여,
컨볼루션 코딩은 현재의 부호어를 만드는데 있어서 현재의 k비트 메시지어만 필요한 것이 아니라 과거의 메시지어도 필요하다. 컨볼루션 코딩은 처리 지연이 작게 되어 실시간 통신이 중요시 되는 경우에 많이 사용된다.

인터리빙(interleaving)이 지향하는 효과는 연속된 오류를 분산시켜서 한 부호어 내에는 오류정정 능력 이하의 비트 오류가 발생되도록 하는 것이다.

채널용량, 정보용량
엔트로피 H(x) 채널입력이 관찰되기 전에 채널입력의 불확실성을 나타내고, 조건부 엔트로피 H(x|y)는 채널출력이 관찰된 후에 채널 불확실성을 표현하므로 H(x)-H(x|y)의 차이는 채널출력을 관찰함으로써 분석되는 채널입력에 대한 불확실성을 표현한다.
그리고 이것을 상호정보량 I(x;y)이라고 한다.


I(x;y)=H(x)-H(x|y)
결합엔트로피 I(x;y)=H(x)+H(y)-H(x,y) , H(x,y)=∑∑p(x,y)*log2(1/(p(x,y))

이산 무기억 채널의 채널용량을 채널의 단일 사용에 있어 최대 평균 상호정보로 정의한다.
C=max[I(x;y)]

채널용량 C는 채널사용당 비트, 혹은 전송당 비트로서 측정된다.
 




2진대칭채널(BSC)

채널은 두 개의 입력심벌과 두 개의 출력심벌을 갖는다. 채널은 0이 보내졌을 때 1로 받을 확률은 1이 보내졌을 때 0으로 받을 확률과 같기 때문에 대칭적이다. 오류의 조건부 확률은 p로 정의되어진다.

H(x)는 채널입력 확률이 1/2일 때, 최대(1)로 된다.

C=1+p*log2(p)+(1-p)*log2(1-p) =1-H(p)

① 채널은 잡음이 없을 때, p=0이라 하면 채널용량 C는 채널사용당 한 비트의 최대값이 되는데 그것은 각 채널입력에서의 정확한 정보가 된다. 위의 p값에서 엔트로피 함수 H(p)는 최소값 0에 도달한다.
② 잡음에 의한 조건부 오류확률 p=1/2일 때 채널용량 C는 최소치 0에 도달하고 반면에 엔트로피 함수 H(p)는 최대값 1에 도달하는데 그러한 경우의 채널은 쓸모 없다.

p=linspace(0,1,100); C=1+p.*log2(p)+(1-p).*log2(1-p); plot(p,C)
 





Shannon's 2st theorem : 2진 대칭채널에서 채널용량 C와 같거나 이하인 어떤 부호율 r(code rate)에 대해서도 원하는 평균 오류확률보다 작은 오류확률을 가지는 부호가 존재하한다.

(예시) 반복부호에서 코드 rate에 대한 평균 오류확률
r=1일때, Pe는 1/100 , r=1/3일때, Pe=∑i=2~3 (3C2)*p^i*(1-p)^(3-i) ∴ Pe=3/10000



가우시안 확률분포함수(pdf) f_X (x)=1/(σ√2π) e^([-(x-m)^2/2σ^2 ])

평균 0, 분산 1로 정규화된 가우시안 누적분포함수(=Q함수)

Q(x)=∫_x~∞[1/√2π e^([-u^2/2]) du]

erfc(x)=2Q(√2x) ※ Q(0)=0.5 , Q(1)=0.1587 , erfc(x)=1-erf(x)

비트오류확률 Pb=Q(sqrt(Eb((1-p)/N0)) , p는 상관계수 ∴ BPSK의 Pb=Q(sqrt(SNR))

Xk가 전송된 신호의 표본이라고 한다. 채널의 출력은 평균이 영이고 전력스펙트럼밀도가 N0/2의 부가백색 가우시안 잡음(AWGN)에 의해 방해를 받는다. 잡음도 B Hz로 대역제한된다. 연속랜덤변수 Yk=Xk+Nk 가 수신된 신호의 표본이라고 하자

잡음표본 Nk는 평균이 영이고, 분산이 σ^2=N0*B , E[Xk^2]=P

C=max[I(Xk;Yk):E[Xk^2]=P] =H(Yk)-H(Yk|Xk)

C=H(Yk)-H(Nk) =1/2*log2(2πe(P+σ^2)-1/2*log2(2πeσ^2) =1/2*log2(1+P/(σ^2))

단위시간당 정보용량 C=B*log2(1+P/(N0*B))
 


Shannon's 3st theorem : 시스템이 전력제한되고, 대역제한된 가우시안 채널에 대해 오류 없는 전송의 속도에 대한 기본적인 제한을 정의
C=B*log2(1+P/(N0*B)) bit/sec

(예시) 대역폭이 19.2[kHz]인 RF채널의 용량이 9.6[kbps]가 되도록 하기 위해 필요한 신호의 SNR은? 9600=19200*log2(1+SNR)[bps]  ∴ SNR=-3.83[dB]


Linear Block Code
Hamming 거리 d는 임의의 두 부호어 벡터의 같은 위치에 있는 원소들 중 서로 다른 원소의 개수로 정의된다.
오류검출능력 q=d_min-1 , 오류정정능력 t=(d_min-1)/2

블록코딩은 k 비트의 message vector를 n 비트의 code vector로 매핑시키는 것을 말하는데, 길이 n의 codeword 2^k개로 구성된 code 를 (n, k) block code라고 한다.
c=mG , G_k*n=[I_k*k|P_k*(n-k)] , H=[P'_(n-k)*k|I_(n-k)*(n-k)] , G*H'=P+P=0
syndrome vector s=r*H'=(c+e)H'=e*H'
c=r-e
디코더는 s에 해당하는 에러유형 e를 찾아서 그것을 수신벡터 r로부터 뺌으로써 원래의 올바른 code vector c를 구한다.




% (7, 4) 블록코딩

SNRbdB=10; MaxIter=1e6;

n=3; N=2^n-1; K=2^n-1-n; % Codeword (Block) length and Message size

Rc=K/N; % Code rate

SNRb=10.^(SNRbdB/10);

SNRbc=SNRb*Rc; % SNRbc*n=SNR*k

sqrtSNRbc=sqrt(SNRbc);

pemb_uncoded=Q(sqrt(SNRb)); % Uncoded msg BER with BPSK

et=Q(sqrt(SNRbc)); % Crossover probability

L=2^K;

for i=1:L
M(i,:)=deci2bin1(i-1,K); % All message vectors

end

P=[1 1 0; 0 1 1; 1 1 1; 1 0 1];

G=[P eye(K)];, H=[P' eye(N-K)]; % [H,G]=Hammgen(n);

Hamming_code=rem(M*G,2);

Min_distance =min(sum(Hamming_code(2:L,:)')); % weight

No_of_correctable_error_bits=floor((Min_distance-1)/2); % 오류정정비트수

E= eye(N); % Error patterns 에러유형

S= rem(E*H',2); NS=size(S,1); % The syndrome matrix

nombe=0;

for iter=1:MaxIter
msg=randint(1,K); % Message vector
coded= rem(msg*G,2); % Coded vector
modulated= 2*coded-1; % BPSK-modulated vector
r= modulated +randn(1,N)/sqrtSNRbc; % Received vector with noise
r_sliced= r>0; % Sliced

r_c=r_sliced; % To be corrected only if it has a syndrome

s= rem(r_sliced*H',2); % Syndrome

for m=1:NS % Error correction depending on the syndrome

if s==S(m,:), r_c=rem(r_sliced+E(m,:),2); break; end

end
notbe1=sum(coded~=r_sliced); notbec1=sum(coded~=r_c);

% 2개이상 비트 오차가 발생했을 때 정정은 오히려 오차개수를 증가시킴

nombe=nombe+sum(msg~=r_c(N-K+1:N)); if nombe>100, break; end

end

pemb=nombe/(K*iter); % Message bit error probability

pemb_t=prob_err_msg_bit(et,N,No_of_correctable_error_bits);

fprintf('\n Uncoded Messsage BER=%8.6f',pemb_uncoded)

fprintf('\nMessage BER=%5.4f(Theoretical BER=%5.4f)\n',pemb,pemb_t)

G, Hamming_code, H', [E S]

 

G =
     1     1     0     1     0     0     0
     0     1     1     0     1     0     0
     1     1     1     0     0     1     0
     1     0     1     0     0     0     1

Hamming_code =
     0     0     0     0     0     0     0
     1     0     1     0     0     0     1
     1     1     1     0     0     1     0
     0     1     0     0     0     1     1
     0     1     1     0     1     0     0
     1     1     0     0     1     0     1
     1     0     0     0     1     1     0
     0     0     1     0     1     1     1
     1     1     0     1     0     0     0
     0     1     1     1     0     0     1
     0     0     1     1     0     1     0
     1     0     0     1     0     1     1
     1     0     1     1     1     0     0
     0     0     0     1     1     0     1
     0     1     0     1     1     1     0
     1     1     1     1     1     1     1

H' =
     1     1     0
     0     1     1
     1     1     1
     1     0     1
     1     0     0
     0     1     0
     0     0     1

[E | S] =
     1     0     0     0     0     0     0     1     1     0
     0     1     0     0     0     0     0     0     1     1
     0     0     1     0     0     0     0     1     1     1
     0     0     0     1     0     0     0     1     0     1
     0     0     0     0     1     0     0     1     0     0
     0     0     0     0     0     1     0     0     1     0
     0     0     0     0     0     0     1     0     0     1


 참고로 BER은 berawgn 과 비교해 보면 큰 향상이 있음.


% BER of a block code

n = 23; k = 12; % Lengths of codewords and messages
dmin = 7; % Minimum distance
EbNo = 1:10;
ber_block = bercoding(EbNo,'block','hard',n,k,dmin);
berfit(EbNo,ber_block) % Plot BER points and fitted curve.
ylabel('Bit Error Probability');
title('BER Upper Bound vs. Eb/No, with Best Curve Fit');




Convolutional Code와 Vitebi Decodeing
컨볼루션 코딩은 정보 비트열이 직렬로 채널 부호화기에 입력되고 전송량이 증가한 부어열이 직렬로 출력되어 전송되는 방식이다. (n, k , K)코드로 표현하는데, K는 n비트의 부호어를 생성하기 위하여 필요한 메시지어의 개수를 나타내고, 이를 구속장이라 한다.

(2, 1, 3) 컨볼루션 부호화기를 예로 살펴보자!
부호화기에 연결벡터는 g1=(111), g2=(101)





trel = poly2trellis(3, [7 5]) % constraint length(구속장) K=3 , 1 1 1 ->7 , 1 0 1 ->5

msg = [1 1 0 1 1 0 0]

[code state]=convenc(msg,trel), tblen = 1; % traceback length

out = vitdec(code, trel, tblen, 'cont', 'hard')


code =

1 1 0 1 0 1 0 0 0 1 0 1 1 1


state =
0

out =
0 1 1 0 1 1 0

비터비 복호 알고리즘은 최우복호를 위한 연산량을 줄일 수 있는 효율적인 방법을 제시한다. 비터비 알고리즘은 각 순간마다 각 상태로 들어오는 트렐리스 경로들과 수신 코드열의 거리를 비교하여 최우경로의 가능성이 없는 경로는 고려 대상에서 제거함으로써 연산량을 줄일 수 있다.

이상,,, 상당히 이해하기 난해한 부분이 많이 있지만, 가능한 감이라도 오도록 참고내용을 요약해서 그대로 적는다. 디지털 데이터전송을 잘하기 위해 비트열을 가지고 이렇게 줄였다 늘렸다 곱했다 나눴다 찌지고 뽂는 가공(프로세싱)을 하는 것을 보면, 정보에는 공장에서 찍어만드는 공산품보다도 더 눈물겨운 수고로움이 축척되어 있음을 알게된다.
마지막으로 비터비란 통신의 대가는 퀄컴의 창업주라고 한다. 역시 자신의 알고리즘을 CDMA 휴대폰에 전송에러를 획기적으로 줄이데 활용했으리라. 뿐만아니라 여러 공학적인 응용에서 비터비알고리즘이 많이 쓰이는 듯 하다. 책에서 보던 공학자가 세계의 몇백번째 부자로 사회적으로 성공하고, 아직도 실존하고 있다는 것이 나에게는 약간은 충격적이고, 부럽다.ㅎ

2011년 11월 7일 월요일

소스코딩

소스코딩은 어떤규칙에 의해서 정보를 압축시켜 전송량을 줄이고자하는 절차이다.

엔트로피는 심벌당 평균 정보량 이다. (불확실성의 측정량)

H(x)=∑p∙log_2 (1/p)

엔트로피는 사건의 확률이 같을때, 최대가 된다.

 (예시)  00, 01, 10, 11 네개의 심볼이 p=1/4 확률로 발생하면, 엔트로피는 2 bits

Shannon's 1st theorem : 정보를 코딩하는데 필요한 bit의 개수는 적어도 그 엔트로피 이상 필요



% 허프만 코딩

p=[0.2 0.15 0.13 0.12 0.1 0.09 0.08 0.07 0.06];    % 어떤 사건이 일어날 확률들

h={'1' '0'};

M=length(p);  N=M-1; p=p(:); % make p a column vector

if M>2

  pp(:,1)=p;

  for n=1:N

     [pp(1:M-n+1,n),o(1:M-n+1,n)]=sort(pp(1:M-n+1,n),1,'descend'); % in descending order

     if n==1, ord0=o; end  % Original descending order

     if M-n>1

       pp(1:M-n,n+1)=[pp(1:M-1-n,n); sum(pp(M-n:M-n+1,n))];

     end

  end

  for n=N:-1:2

     tmp=N-n+2; oi=o(1:tmp,n);

     for i=1:tmp, ht{oi(i)}=h{i}; end

     h=ht; h{tmp+1}=h{tmp};

     h{tmp}=[h{tmp} '1'];

     h{tmp+1}=[h{tmp+1} '0'];

  end

  for i=1:length(ord0), ht{ord0(i)}=h{i}; end

  h=ht;

end

L=0;

for n=1:M

  L=L+p(n)*length(h{n});% Average codeword length

end

H=-sum(p.*log2(p)); % Entropy


>> L,H,pp,o,h
L =    3.1000
H =     3.0731

pp =
0.2000 0.2000 0.2000 0.2200 0.2600 0.3200 0.4200 0.5800
0.1500 0.1500 0.1700 0.2000 0.2200 0.2600 0.3200 0.4200
0.1300 0.1300 0.1500 0.1700 0.2000 0.2200 0.2600 0
0.1200 0.1300 0.1300 0.1500 0.1700 0.2000 0 0
0.1000 0.1200 0.1300 0.1300 0.1500 0 0 0
0.0900 0.1000 0.1200 0.1300 0 0 0 0
0.0800 0.0900 0.1000 0 0 0 0 0
0.0700 0.0800 0 0 0 0 0 0
0.0600 0 0 0 0 0 0 0

o =
1 1 1 6 5 4 3 2
2 2 7 1 1 1 1 1
3 3 2 2 2 2 2 0
4 8 3 3 3 3 0 0
5 4 4 4 4 0 0 0
6 5 5 5 0 0 0 0
7 6 6 0 0 0 0 0
8 7 0 0 0 0 0 0
9 0 0 0 0 0 0 0

h =
Columns 1 through 8
'00' '110' '101' '011' '010' '1111' '1110' '1001'
Column 9
'1000'

부연설명을 하자면, 1 ~9까지의 심벌의 발생확률이 p라고 했을때,
100심벌을 전송한다고 가정해 보자!
메시지를 2진코드로 4bit씩 그대로 코딩하면, 400bit 가 필요한 반면,
허프만코딩은 20*2+(15+13+12+10)*3+(9+8+7+6)*4=310bit 로 정보를 압축할 수 있다.

이와달리, 심벌들의 발생과 확률분포가 필요없이 주어진 소스 시퀀스에 대한 codetable을 만들어 사용하는 Lempel-Ziv-Welch 코딩도 있다.
Lempel-Ziv 알고리즘은 허프만을 대신해 파일압축을 위한 기본으로 널리 이용되고 있다.

지금까지는 통계적 중복성 제거 를 통한 정보압축이었다.


시각적/청각적 특징 을 기반한 정보의 압축으로는
영상의 디지털화 형식(원래 RGB=8:8:8 → Ycbcr=8:2:2, 8:4:4), 음성의 마스킹 효과가 있다.


시간적 중복성 제거
영상에서 움직임 예측 및 움직임 보상 이론


공간적 중복성 제거
예로는, JPEG의 압축기술인 DCT를 들 수 있다.

DCT(Discrete Cosine Transform)은 직교 변환을 이용한 것으로, 공간상의 영상 데이터를 주파수상의 에너지 분포로 변환시키는 것이다.

F(u,v)=(1/4)∙C(u)C(v)∑_(i=0~7)∑_(j=0~7)f(i,j)cos((2i+1)uπ/16) cos((2j+1)vπ/16)
여기서 C(x)={(1/√2 ,x=0 @ 1 ,otherwise)}



Xr=[]; X=[];
for u=0:7
           for v=0:7
                     for i=0:7
                                for j=0:7
                    x(i+1,j+1)=cos((2*i+1)*u*pi/16).*cos((2*j+1)*v*pi/16);
                                end
                     end
                     Xr=[Xr x];
                subplot(8,8,u*8+v+1), imshow(x,[-1,1])
           end
           X=[X; Xr]; Xr=[];
end
 



64개의 DCT 기본함수


용량큰 JPEG를 화면에 뿌려줄때, 흐릿하다가 좀 뒤에 또렷해지는 것은...
위에 저주파의 신호부터 지그재그로 읽어들이는 도중에 먼저 화면에 뿌리고, 고주파의 신호까지 완벽히 읽은 뒤에 마지막으로 다시 화면에 뿌려주기 때문일 것이다.

부연설명을 하자면, 영상을 DCT변환하면 수직·수평 저주파 계수들의 값은 크고, 고주파로 갈 수록 값이 작아지는데, 이를 양자화하는데서 압축을 하면 정보량을 많이 줄일 수 있다.
 

RGB = imread('autumn.tif');
        I = rgb2gray(RGB);
        J = dct2(I);
        imshow(log(abs(J)),[]), colormap(jet), colorbar

        J(abs(J)<10) = 0;
        K = idct2(J);
        figure, imshow(I)

        figure, imshow(K,[0 255])
 
 


웨이블릿 변환

 정지영상압축 표준인 JPEG2000에서 매우 중요한 부분을 차지하고 있으며, MPEG-4, H.264등에도 사용되고 있다. 또한 잡음제거, 에지 검출과 같은 다양한 영상신호 처리에 응용된다.
웨이블릿 변환에 사용되는 기저 함수의 집합은 하나의 기본 웨이블릿 기저 함수에 대한 시간축 방향으로의 확대 및 축소 그리고 평행 이동을 통해 얻어진다. 기본 웨이블릿 기저 함수는 특별한 형태의 band-pass 필터로 생각할 수 있다. 웨이블릿 변환에서는 주파수 대역이라는 용어 대신에 스케일(scale)이라는 용어를 주로 사용한다.


Ψ_j,k(t) : 웨이블릿 크기 압축계수(scale, j), 시간축으로의 이동 전이계수(translattion, k)

DWT : f(t)=∑_(j⊂Z)C_0*Φ_0,j(t)+∑_(k≥0))∑_(j⊂Z)d_k,j*Ψ_k,j(t)



load wbarb;
wavelet = 'haar'; % Define wavelet of your choice
level = 2; % Define wavelet decomposition level
[C S] = wavedec2(X,level,wavelet); % Compute multilevel 2D wavelet decomposition
% decomposition vector C = [ A(N) | H(N) | V(N) | D(N) | ... | H(1) | V(1) | D(1) ]
% the corresponding bookkeeping matrix S


% Define colormap and set rescale value
colormap(map); rv = length(map);
% Plot wavelet decomposition using square mode

A = cell(1,level); H = A; V = A; D = A;
for k = 1:level
A{k} = appcoef2(C,S,wavelet,k); % Extract 2-D approximation coefficients
[H{k} V{k} D{k}] = detcoef2('a',C,S,k); % Extract 2-D detail coefficients
A{k} = wcodemat(A{k},rv); % Extended pseudocolor matrix scaling % approximation coefficients
H{k} = wcodemat(H{k},rv); % hori. detail coefficients
V{k} = wcodemat(V{k},rv); % vert. detail coefficients
D{k} = wcodemat(D{k},rv); % diag. detail coefficients
end
dec = cell(1,level);
dec{level} = [A{level} H{level} ; V{level} D{level}];

for k = level-1:-1:1
dec{k} = [imresize(dec{k+1},size(H{k})) H{k} ; V{k} D{k}];
end
figure(1);
image(dec{1});
title(['Decomposition at level ',num2str(level)]);

% Peak Signal To Noise Rate
X0 = waverec2(C,S, wavelet);
error=X-X0;
s_error=error^2;
ms_error=sum(sum(s_error),2)/(256*256);
Peak_SNR=10*log(255*255/ms_error)

% Analyze Frequency Characteristic of above wavelet filters
[D_L, D_H, R_L, R_H] = wfilters(wavelet);
den = [1];
figure(2);
freqz(D_L, den);
figure(3);
freqz(D_H, den);
figure(4);
subplot(221); stem(D_L);
title('Decomposition low-pass filter');
subplot(222); stem(D_H);
title('Decomposition high-pass filter');
subplot(223); stem(R_L);
title('Reconstruction low-pass filter');
subplot(224); stem(R_H);
title('Reconstruction high-pass filter');
xlabel(['The four filters for ' wavelet])
 

 




마찬가지로,,, 대부분의 정보는 왼쪽 위에 모서리에 집중되어 있다.

참고로, wfilters의 종류에 따른 웨이블릿 필터의 주파수응답(Fig2, 3)과 임펄스응답(Fig4)을 확인 할 수도 있다.

2011년 10월 1일 토요일

CATV (케이블방송)

현재 IPTV와 위성방송을 위시한 여러 신기술과의 치열한 경쟁속에서도 케이블TV는 꿋꿋이 잘 견디며, 세력을 확장하는 것 같다. 개인적으로 케이블TV는 서민의 느낌으로 더 친밀감이 든다. 또 요새 다소 선정적이지만 재미있는 PP(Program Provider)프로그램도 많다. '막돼먹은 영애씨', '슈스케', '롤로코스터' 등을 보면, 지상파보다도 더 케이블이 대세이다.
동네 허접한 장비 갖다놓고 비디오를 틀어 주던 중계유선방송사에서 시작해서 몸집을 키우면서 지금의 종합유선방송사(SO, System Operator)까지, 케이블방송을 보면, 국내 열악한 방송통신계의 흥망성쇠의 긴 역사가 고스란이 배여 있다고 생각한다. 케이블방송은 저렴하면서 인터넷, 방송, VOD, 게임, 인터넷전화 등 다양한 서비스가 가능하다. 우리는 가까이 있으면서 편하게 즐기는 케이블방송을 무관심해서인지, 잘 모르는 것 같다. 오늘은 알면 더 재밌는 케이블방송에 대해서 알아 보자~!!


케이블방송은 행정구역 상 '구' 정도 급에 한개 업체가 독점하거나 대부분 두개 이상 업체가 경쟁을 하고 있다. 대기업(CJ, 현대백화점, 티브로드, C&M)들이 독립적인 SO들을 전국적으로 여러게 보유하고 있으며, 독과점은 법(전체 가입자 1/3 이하)으로 제한되어 있다.

전체적이고 개괄적인 구성은 이렇다. (광 케이블이랑 동축케이블의 혼합)


CATV 방송망 Headend 다이어그램

먼저 스튜디오에서 제작한 프로그램을 PP 나 DMC(Digital Media Center)에서 제공받고, KT나 전국망사업자들에 의해 광케이블이나 위성으로 케이블방송국의 Headend 까지 공급받는다.
Headend에서는 몇개의 중심 분배센터를 거쳐 개별 셀까지 광신호를 전송한다. 이후부터는 동축케이블이 이용되는데, 우리가 한전전주에서 송전선 아래서 보는 것이 그것이다. 물론 인터넷서비스공급업체의 망과도 짬봉되기도해서 정리안되고 지저분해 보일 수 있다.
이제 전주에 걸린 TAB단에서 집으로 역시 동축케이블만 연결을 하면, 드디어 TV도 보고, 인터넷도 하고, 인터넷전화도 쓸 수 있다!
여기서, CMTS는 서비스가입자를 확인하고 기록하며, 대역폭을 할당해서 인터넷전송속도를 제어하는 역할을 한다.


아날로그방송


사실 2012, 내년 이면, 지상파에서 없어질 기술이다. 케이블은 계속 서비스 할 것이지만,,,
그래도 자세하게 잘 설명된 링크를 걸어 둔다.

http://ko.wikipedia.org/wiki/NTSC

국내에서 아날로그방송의 표준은 NTSC이고, 1ch 대역폭이 6MHz씩, 1GHz이하 주파수를 이용하고 있다. 영상은 VSB 변조, 음성은 FM 변조된다.
대략 디지털 정보량으로 환산하면, 전송률(오디오 제외)은,,,

525(주사선)*858(주사선 당 화소수)*30(60Hz 비월주사)*16(Y=CbCr=8비트)=216Mbps

아날로그TV 영상을 뿌려주기 위해서, 위 전송률을 텔레비전이 매초 처리해야 된다.


컴포넌트(RGB), NTSC 컴포지트(YCrCb) 신호 파형


RGB = imread('naver.bmp');
Ycbcr = rgb2ycbcr(RGB); Y=Ycbcr(:,:,1);
subplot(2,1,1), image(uint8(RGB)), title('RGB')
subplot(2,1,2), imshow(uint8(Y)), title('luminance Y')

 

matlab 에서 RGB를 Ycbcr로 변환

중요한것은, NTSC는 기존 흑백TV와의 호환을 우선한 기술이다.
또한 휘도 신호(luminance, Y)와 색차 신호(chrominance, 임의의 색과 그 색과 같은 휘도를 가진 기준색과의 색적인 차, Cb,Cr->I,Q)가 함께 있는 컴포지트 인터페이스이다. 색차 신호는 컬러 부반송파에 의해 휘도 신호 반송파 보다 높은 주파수에 실리게 된다.
컬러부호화에서 부반송파에 실린 I 와 Q 신호의 위상(phase) 관계는 TV 카메라에 잡힌 색상의 색도(color hue)와 관련이 있고, 이들의 진폭(amplitude)은 원 신호의 색상의 채도(saturation, purity)와 관련이 있다고 한다.


NTSC 컴포지트(YCrCb) 신호 파형의 스펙트럼


%--------------------------------------------------------------------------------------------------
%          simulates FM for the message signal
%--------------------------------------------------------------------------------------------------
Ac=1; fc=50; wc=2*pi*fc;              % Amplitude/Frequency of carrier
Tb=0.1;                                            % Bit interval time
T=1/fc/8; Fs=1/T;                            % Sampling period, frequency
Nb=Tb/T; n=2^(nextpow2(3*Nb)); t=[1:n]*T;         
                                                        % Nb=Tb 내 샘플수, n=총 샘플수, time vector
f =[-Fs/2: Fs/n: Fs/2];                      % Freq. vector

m= ones(Nb,1)*[4 -8 -4]; m=m(:).';                    % message signal m(t)
m=[m, zeros(1,n-length(m))]; x=m;
Mf=fftshift(fft(m)); Mf=[Mf Mf(1)]*T;              % Spectrum of Message signal m(t)

dm(1)=0; T2=T/2;
for i=1:n-1
 dm(i+1)=dm(i)+(m(i)+m(i+1))*T2;
end
kf=30;                                                            % deviation constant , th=kf*m(t)*T
s=Ac*cos(wc*t+kf*dm);
Sf=fftshift(fft(s)); Sf=[Sf Sf(1)]*T;              % Spectrum of modulated signal s(t)
s_th=angle(Sf);
sa=hilbert(s);                                    
% xa=hilbert(x) , xa=x+i*x_hat=x_l*exp(j*w1*t)=xc+i*xs <-> Xaf(w)=2*u(w)*Xf(w) 
% abs(xa)= x envelope
ssb=real(hilbert(x)).*cos(wc*t)-imag(hilbert(x)).*sin(wc*t); SSB_f=abs(fftshift(fft(ssb)));
Sa=fftshift(fft(sa)); Sa=[Sa Sa(1)]*T;
sa_th=angle(sa);
th=unwrap(angle(hilbert(s)))-wc*t;   % phase of analytic signal
% th=tan-1(xs/xc)=tan-1(x_hat/x)-wt ,  
% x=xc*cos(wt)-xs*sin(wt) , x_hat=xs*cos(wt)+xc*sin(wt)
th./dm

y=[0 diff(th)/T/kf];
Yf=fftshift(fft(y)); Yf=[Yf Yf(1)]*T;              % Spectrum of modulated signal s(t)

subplot(4,3,1), plot(t, m), title('Message signal m(t)')
subplot(4,3,2), plot(f, abs(Mf)), title('Spectrum of message')
subplot(4,3,4), plot(t, s), title('FM modulated signal')
subplot(4,3,5), plot(f, abs(Sf)), title('Spectrum of modulated signal')
subplot(4,3,6), plot(t, unwrap(angle(hilbert(s))), t, wc*t, 'r'), title('phase of FM modulated signal')
subplot(4,3,7), plot(t, abs(sa)), title('Envelope of FM modulated signal') 
subplot(4,3,8), plot(f, abs(Sa)), title('Spectrum of analytic signal')
subplot(4,3,9), plot(t, sa_th), title('phase of analytic signal')
subplot(4,3,10), plot(t, y), title('Demodulated signal')
subplot(4,3,11), plot(f, abs(Yf)), title('Spectrum of Demodulated signal')





디지털방송

국내에서 디지털방송의 표준은 ATSC이고, full HD는 주사선은 1080i 이다.
지상파 DTV는  6MHz 대역폭으로 , MPEG 2 인코딩하고 8VSB 변조되어  TS bit rate 19.393Mbps 이다. 디코딩하면 훨씬 정보량이 많을테지만...
(참고로 인코딩 및 다중화 방식 중에서, 저장매체에 저장하기 위한 프로그램 스트림 PS-program stream 이 있고, 네트워크에서 전송 또는 디지털 방송을 위한 트랜스포트 스트림 TS-transport stream이 있다.)
케이블방송은 옥상에 설치된 V/UHF안테를 통하여 지상파DTV를 직접수신하여, 가입자들에게 동축케이블을 통해 그대로 8VSB 전송한다.  물론 재생해서 전력을 확 줄여 보내겠지만...

ATSC 디지털방송의 8VSB 스펙트럼
또 지상파DTV 재송신과는 다르게, 다른 채널의 디지털 방송은 대부분 MPEG 4 인코딩하고 256 QAM변조해서 송출된다. 여기서 부터는 고화질이라서 유료서비스이다.

256 QAM 변조한 한 채널의 스펙트럼

디지털방송은 아날로그보다 화질과 음성이 훨씬 좋은데, 한 채널의 대역폭(6MHz)은 동일한 것을 보면, 그만큼 부호화(MPEG 인코딩)하는 기술이 놀랍고, 그 혜택이 크다고 볼 수 있다.


인터넷


케이블방송은 인터넷서비스도 제공하는데, 케이블방송의 상향채널을 통해 가입자로 부터 데이터를 올려받고, 하향채널을 통해 가입자들에게 데이터를 전송해 준다. 핵심적인 장비로 CMTS가 있는데, 인터넷에서 발췌했다.


CMTS (cable modem termination system)는 케이블 모뎀(Cable Modem) 데이터를 인터넷으로 전송하기 위한 데이터용 패킷으로 바꾸어주는 장비이다. 
CMTS는 데이터 네트워크에서 RF케이블 네트워크로 라우팅 기능을 제공하는 CAR(Cable Access Router)카드와 하향 데이타를 위한 64/256QAM을 처리할 수 있는 QAM Modulator, 상향 데이타를 위한 QPSK/QAM Modulator로 구성 가능하다.
케이블모뎀종단치(CMTS) HFC 망이 끝나는 분배센터 내에 위치하여 양방향 HFC망을 통해 케이블 모뎀장비와 인터페이스를 수행하여 외부망과 연결된다. CMTS는 케이블 모뎀과 인터페이스를 위해 다양한 종류의 케이블 모뎀 카드 지원 및 케이블 모뎀 인증을 담당하며, 상 하양 채널의 주파수를 지정하고 채널에 대한 데이터를 암호호한다. CMTS는 대게 라우터, 변조기, 라인카드가 조합된 형태로 구성이 되며, 인터넷/인트라넷 백본 데이터 네트워크와 지역 액세스 케이블 네트워크 간에 고속 통신을 가능하게 한다. 일부 벤더들은 L#라우터 대신에 L2 브리징과 스위칭 기술을 사용하여 CMTS와 분산된 라우터를 연결한다.


케이블방송 인터넷서비스 다이어그램(6MHz RF 주파수 4ch 하향신호)


막상 정리하고 보니, 인터넷에서 너무 퍼오고, 짜집기가 난무한거 같다. 저작권도 약간은 걱정되기도 한다. 또 복잡한 주제에 대해 처음에는 거창하게 시작했지만, 내공의 부족함을 느끼고 대략 마무리 한다.^^
그래도 나름의 의미를 찾자면, 깊이는 부족하지만 폭 넓게 한눈에 겉 핥기식으로 볼수 있다는 점이다.ㅎㅎㅎ
차차 좀 더 자료는 보충할 생각이다~! 이상~