2010年12月15日 星期三

PCA eigenface

 出處:http://www.pages.drexel.edu/~sis26/Eigencode.htm

% read and show image
M=4;

% Chosen std and mean.
% It can be any number that it is close to the std and mean of most of the images.
um=100;
ustd=80;

S=[];    % img matrix
figure(1);
for i=1:M
    str=strcat(int2str(i),'fr','.bmp');    % concatenates two strings that form the name of the image
    eval('img=rgb2gray(imread(str));');
    subplot(ceil(sqrt(M)),ceil(sqrt(M)),i)
    imshow(img)
    if i==3
        title('Training set','fontsize',18)
    end
drawnow;
[irow icol]=size(img);    % get the number of rows (N1) and columns (N2)
temp=reshape(img',irow*icol,1);    % creates a (N1*N2)x1 vector
S=[S temp];    % S is a N1*N2xM matrix after finishing the sequence
end

% Here we change the mean and std of all images. We normalize all images.
% This is done to reduce the error due to lighting conditions and background.
for i=1:size(S,2)
    temp=double(S(:,i));
    m=mean(temp);
    st=std(temp);
    S(:,i)=(temp-m)*ustd/st+um;
end

% show normalized images
figure(2);
for i=1:M
    str=strcat(int2str(i),'.jpg');
    img=reshape(S(:,i),icol,irow);
    img=img';
    eval('imwrite(img,str)');
    subplot(ceil(sqrt(M)),ceil(sqrt(M)),i)
    imshow(img)
    drawnow;
    if i==3
        title('Normalized Training Set','fontsize',18)
    end
end

% mean image
m=mean(S,2);  % obtains the mean of each row instead of each column
tmimg=uint8(m); % converts to unsigned 8-bit integer. Values range from 0 to 255
img=reshape(tmimg,icol,irow); % takes the N1*N2x1 vector and creates a N1xN2 matrix
img=img';
figure(3);
imshow(img);
title('Mean Image','fontsize',18)


% Change image for manipulation
dbx=[];    % A matrix
for i=1:M
    temp=double(S(:,i));
    dbx=[dbx temp];
end

%Covariance matrix C=A'A, L=AA'
A=dbx';
L=A*A';
% vv are the eigenvector for L
% dd are the eigenvalue for both L=dbx'*dbx and C=dbx*dbx';
[vv dd]=eig(L);
% Sort and eliminate those whose eigenvalue is zero
v=[];
d=[];
for i=1:size(vv,2)
    if(dd(i,i)>1e-4)
        v=[v vv(:,i)];
        d=[d dd(i,i)];
    end
end

%sort, will return an ascending sequence
[B index]=sort(d);
ind=zeros(size(index));
dtemp=zeros(size(index));
vtemp=zeros(size(v));
len=length(index);
for i=1:len
    dtemp(i)=B(len+1-i);
    ind(i)=len+1-index(i);
    vtemp(:,ind(i))=v(:,i);
end
d=dtemp;
v=vtemp;

%Normalization of eigenvectors
for i=1:size(v,2) %access each column
    kk=v(:,i);
    temp=sqrt(sum(kk.^2));
    v(:,i)=v(:,i)./temp;
end

%Eigenvectors of C matrix
u=[];
for i=1:size(v,2)
    temp=sqrt(d(i));
    u=[u (dbx*v(:,i))./temp];
end

%Normalization of eigenvectors
for i=1:size(u,2)
    kk=u(:,i);
    temp=sqrt(sum(kk.^2));
    u(:,i)=u(:,i)./temp;
end

% show eigenfaces
figure(4);
for i=1:size(u,2)
    img=reshape(u(:,i),icol,irow);
    img=img';
    img=histeq(img,255);
    subplot(ceil(sqrt(M)),ceil(sqrt(M)),i)
    imshow(img)
    drawnow;
    if i==3
        title('Eigenfaces','fontsize',18)
    end
end

% Find the weight of each face in the training set
omega = [];
for h=1:size(dbx,2)
    WW=[];
    for i=1:size(u,2)
        t = u(:,i)';
        WeightOfImage = dot(t,dbx(:,h)');
        WW = [WW; WeightOfImage];
    end
    omega = [omega WW];
end

% Acquire new image
% Note: the input image must have a bmp or jpg extension.
% It should have the same size as the ones in your training set.
% It should be placed on your desktop
InputImage = input('Please enter the name of the image and its extension \n','s');
InputImage = rgb2gray(imread(InputImage));
figure(5)
subplot(1,2,1)
imshow(InputImage); colormap('gray');title('Input image','fontsize',18)
InImage=reshape(double(InputImage)',irow*icol,1);
temp=InImage;
me=mean(temp);
st=std(temp);
temp=(temp-me)*ustd/st+um;
NormImage = temp;
Difference = temp-m;

p = [];
aa=size(u,2);
for i = 1:aa
    pare = dot(NormImage,u(:,i));
    p = [p; pare];
end
ReshapedImage = m + u(:,1:aa)*p; %m is the mean image, u is the eigenvector
ReshapedImage = reshape(ReshapedImage,icol,irow);
ReshapedImage = ReshapedImage';
%show the reconstructed image.
subplot(1,2,2)
imagesc(ReshapedImage); colormap('gray');
title('Reconstructed image','fontsize',18)

InImWeight = [];
for i=1:size(u,2)
t = u(:,i)';
WeightOfInputImage = dot(t,Difference');
InImWeight = [InImWeight; WeightOfInputImage];
end

ll = 1:M;
figure(68)
subplot(1,2,1)
stem(ll,InImWeight)
title('Weight of Input Face','fontsize',14)

% Find Euclidean distance
e=[];
for i=1:size(omega,2)
q = omega(:,i);
DiffWeight = InImWeight-q;
mag = norm(DiffWeight);
e = [e mag];
end

kk = 1:size(e,2);
subplot(1,2,2)
stem(kk,e)
title('Eucledian distance of input image','fontsize',14)

MaximumValue=max(e)  % maximum eucledian distance
MinimumValue=min(e)    % minimum eucledian distance


2010年8月4日 星期三

使用 AAM-API

先按照doc安裝

安裝完就會有預先編譯好的執行檔(aamc.exe aamcm.exe)

aamc 針對灰階圖
aamcm 針對彩色圖(muti)

開啟aamc.sln重新編譯本來是出現aamcd.exe

不過速度沒有aamc.exe來的快

所以更改編譯模式 使用Release 3 band AAMCM
在重新編譯會比較快
若有error出現重複定義的情況 可能必須取消 連結器->輸入->其他相依性 的 aam-apimDB.lib

之後就可以修改原始碼囉~

2010年5月25日 星期二

Compressive Sensing (matlab code using DCT)

clc;clear
load DCTq;  %讀量化表格
% Read the input image;           
xtatal = double(imread('lena512.bmp')); %讀圖檔 目前應該是黑白
figure(1);
imshow(uint8(xtatal));  %秀出原始圖檔
% keep an original copy of the input signal
for i=0:size(xtatal,1)/8-1                  %每格寬8像素 總共 (寬/8) 格 需整除
    for j=0:size(xtatal,1)/8-1              %每格高8像素 總共 (高/8) 格 需整除
        x0 = xtatal(1+i*8:8+i*8,1+j*8:8+j*8);       %一個block 處理
%         xmean = mean(x0(:));
%         x0=x0-xmean;
        % figure(1);
        % imshow(uint8(x0));
        f=floor((dct2(x0)+floor(Q/2))./Q);      %一個block作DCT和量化後的 係數(稀疏)
        f = f(:);                               %拉成一維 (n*1)
        n= size(x0,1)*size(x0,1);               %總點數 n

        s=length(find(floor((dct2(x0)+floor(Q/2))./Q)~=0)); %不等於0的有s個
       
        M=32;                                           %設定M大小 每個block 都用一樣
%         M=s*log(n/s)+1;                               %需大於s*log(n/s) 公式?
%         if (M<32)                                     %補不足的點(自己設定的)
%             M=32;
%         end

        A = get_A_random(n,M);          %sensing matrix (M*n)
        Psi=dctmtx(n);                  %DCT基底 (n*n)
        y = A*Psi*f;                    %measurment (壓縮後大小 y(M*1)) (編碼)

        %用(l1 magic)tool 求l1 最佳稀疏解 (解碼)
%         % S o l v e u s ing l 1 magic .
%         path(path ,'./Optimization') ;
%         xinit = pinv(A)*y; % i n i t i a l g u e s s = min ene r g y
%         tic
%         xp = l1eq_pd ( xinit , A, [ ] , y , 1e-3);
%         toc

        %用CVX tool 求l1 最佳稀疏解 (解碼)
        cvx_begin
            variable xp(n);
            minimize (norm( xp , 1 ) ) ;
            subject to
                A*Psi*xp == y ;
        cvx_end
       
        %norm(f-xp)/norm(f);        %計算誤差??
        hat_x=reshape(round(xp),8,8);       %(n*1)一維轉換成8*8
        xhatall(1+i*8:8+i*8,1+j*8:8+j*8)=(idct2(floor(hat_x.*Q)));  %反量化&反DCT 還原圖檔
       
        %調整超過的像素值
        [r1,c1]=(find(xhatall>255));
        for l=1:length(r1)
            xhatall(r1(l),c1(l))=255;
        end
        [r2,c2]=(find(xhatall<0));
        for l=1:length(r2)
            xhatall(r2(l),c2(l))=0;
        end

        %不確定對不對的方式
%         xhatall(1+i*8:8+i*8,1+j*8:8+j*8)=reshape(Psi*reshape(hat_x.*Q,64,1),8,8);
    end
end

figure(2);
imshow(uint8(xhatall));     %秀出還原的圖檔

2010年5月24日 星期一

android porting note

參考http://www.mask.org.tw/data/BringUp_Android_on_PXA270.pdf
必須安裝g++ bison flex 等套件

最後遇到
prebuilt/linux-x86/sdl/include/SDL/SDL_syswm.h:75: error: syntax error before '}' token
類似這種問題

所以google大神說
(http://webcache.googleusercontent.com/search?q=cache:1mLuKdvgfvYJ:dev.firnow.com/course/6_system/linux/Linuxjs/20090314/161469.html+prebuilt/linux-x86/sdl/include/SDL/SDL_syswm.h+error:+syntax+error+before+%27*%27+token&cd=1&hl=zh-TW&ct=clnk&gl=tw&client=firefox-a)

以下不確定什麼要裝 (管他的 全裝)
apt-get install sun-java6-jdk
apt-get install flex bison gperf libesd0-dev libwxgtk2.6-dev build-essential zip curl libncurses5-dev zlib1g-dev
apt-get install libsdl-dev

跑快1小時就可以過了

2010年5月18日 星期二

Local Binary Patterns, LBP

Ojala 等人提出的區域二元特徵(Local Binary Patterns, LBP)已在多篇論文中被證實在描述影像上的表現很好[3],LBP 是一種用來描述區域紋理變化的特徵計算方式。

優點是運算簡單、速度也相當的快,且不受陰影困擾,非常適合使用在真實的即時系統(real-time system)上

缺點是在平滑影像(如轉灰階後,灰階值太過於相近的天空及海洋影像)上的效果較差強人意,但既然是描述影像中的紋理資訊,所以也極少會將區域二元特徵應用在平滑影像上。

2010年5月17日 星期一

uniform patterns

uniform patterns,即8-bit 二元編碼中,
最多只能有2 個0 到1 或1 到0 的轉換,例如:00001111、10001111 和11111111
都是uniform patterns

2010年5月11日 星期二

BCB轉灰階code範例

Bitmap *result = new Bitmap();
result->PixelFormat = pf24bit;
result->Width = Image1->Picture->Bitmap->Width;
result->Height = Image1->Picture->Bitmap->Height;
for (int row = 0; row <>Picture->Bitmap->Height ; row++) {
   BYTE *Pixel1=(BYTE *)Image1->Picture->Bitmap->ScanLine[row];
   BYTE *Pixel2=(BYTE *)result->ScanLine[row];
   for (int col = 0; col <>Picture->Bitmap->Width ; col++) {
      Byte Gray = (*Pixel1+*(Pixel1+1)+*(Pixel1+2))/3;
      *Pixel2=*(Pixel2+1)=*(Pixel2+2) = Gray;
      Pixel1+=3;
      Pixel2+=3;
   }
}
Image2->Picture->Bitmap->Assign(result);
delete result;