pmain.m

主要用于图像识别中

tic
% pmain
% get all the models

% mod1=wwx([997,1093,2012,2080]);
% mod2=wwx([1729,1824,2259,2326]);
% mod3=wwx([1289,1386,2750,2819]);
% mod4=wwx([1876,1970,2786,2853]);
% mod5=wwx([2021,2117,3031,3098]);
% mod6=wwx([2314,2409,2541,2606]);
% mod7=wwx([2168,2263,2669,2734]);
% mod8=wwx([2022,2118,1815,1882]);
% mod9=wwx([3190,3288,2128,2209]);
% mod10=wwx([3629,3726,1392,1472]);
% mod11=wwx([3336,3433,2538,2618]);
tic
mod1=wwx([1010,1080,2012,2080]);
mod2=wwx([1740,1811,2259,2326]);
mod3=wwx([1303,1373,2750,2819]);
mod4=wwx([1887,1957,2786,2853]);
mod5=wwx([2033,2104,3031,3098]);
mod6=wwx([2327,2397,2540,2619]);
mod7=wwx([2180,2250,2669,2734]);
mod8=wwx([2035,2106,1802,1882]);
mod9=wwx([3203,3274,2128,2209]);
mod10=wwx([3642,3712,1392,1472]);
mod11=wwx([3349,3420,2538,2618]);
mod12=wwx([1742,1814,2130,2212]);

modstr(1).va=edge(mod1,'canny');
modstr(1).na=1;
modstr(2).va=edge(mod2,'canny');
modstr(2).na=2;
modstr(3).va=edge(mod3,'canny');
modstr(3).na=3;
modstr(4).va=edge(mod4,'canny');
modstr(4).na=4;
modstr(5).va=edge(mod5,'canny');
modstr(5).na=5;
modstr(6).va=edge(mod6,'canny');
modstr(6).na=6;
modstr(7).va=edge(mod7,'canny');
modstr(7).na=7;
modstr(8).va=edge(mod8,'canny');
modstr(8).na=8;
modstr(9).va=edge(mod9,'canny');
modstr(9).na=9;
modstr(10).va=edge(mod10,'canny');
modstr(10).na=10;
modstr(11).va=edge(mod11,'canny');
modstr(11).na=11;
modstr(12).va=edge(mod12,'canny');
modstr(12).na=12;

xx=length(modstr);
sse=strel('square',3);
ccov=zeros(xx,xx);
for ii=1:xx-1
    a=modstr(ii).va;
    [ma,na]=size(a);
    ac=imdilate(a,sse);
    for jj=ii+1:xx
        b=modstr(jj).va;
        [mb,nb]=size(b);
        bc=imdilate(b,sse);
        xmax=0;
        if ma>=mb & na>=nb
            sumb=sum(sum(b));
            for kk=1:ma-mb+1
                for ll=1:na-nb+1
                    outf=and(a(kk:kk+mb-1,ll:ll+nb-1),bc);
                    soutf=sum(sum(outf));
                    forw=soutf/sum(sum(a(kk:kk+mb-1,ll:ll+nb-1)));
                    outb=and(ac(kk:kk+mb-1,ll:ll+nb-1),b);
                    soutb=sum(sum(outb));
                    backw=soutb/sumb;
                    if xmax<min(forw,backw);
                        xmax=min(forw,backw);
                        corcell{ii,jj}=outf;
                        corcell{jj,ii}=outb;
                    end
                end
            end
        elseif ma>=mb & na<nb
            for kk=1:ma-mb+1
                for ll=1:nb-na+1
                    outf=and(a(kk:kk+mb-1,:),bc(:,ll:ll+na-1));
                    soutf=sum(sum(outf));
                    forw=soutf/sum(sum(a(kk:kk+mb-1,:)));
                    outb=and(ac(kk:kk+mb-1,:),b(:,ll:ll+na-1));
                    soutb=sum(sum(outb));
                    backw=soutb/sum(sum(b(:,ll:ll+na-1)));
                    if xmax<min(forw,backw);
                        xmax=min(forw,backw);
                        corcell{ii,jj}=outf;
                        corcell{jj,ii}=outb;
                    end
                end 
            end
        elseif ma<mb & na>=nb
            for kk=1:mb-ma+1
                for ll=1:na-nb+1
                    outf=and(a(:,ll:ll+nb-1),bc(kk:kk+ma-1,:));
                    forw=sum(sum(outf))/sum(sum(a(:,ll:ll+nb-1)));
                    outb=and(ac(:,ll:ll+nb-1),b(kk:kk+ma-1,:));
                    backw=sum(sum(outb))/sum(sum(b(kk:kk+ma-1,:)));
                    if xmax<min(forw,backw);
                        xmax=min(forw,backw);
                        corcell{ii,jj}=outf;
                        corcell{jj,ii}=outb;
                    end
                end 
            end
        elseif ma<mb & na<nb
            suma=sum(sum(a));
            for kk=1:mb-ma+1
                for ll=1:nb-na+1
                    outf=and(a,bc(kk:kk+ma-1,ll:ll+na-1));
                    forw=sum(sum(outf))/suma;
                    outb=and(ac,b(kk:kk+ma-1,ll:ll+na-1));
                    backw=sum(sum(outb))/sum(sum(b(kk:kk+ma-1,ll:ll+na-1)));
                    if xmax<min(forw,backw);
                        xmax=min(forw,backw);
                        corcell{ii,jj}=outf;
                        corcell{jj,ii}=outb;
                    end
                end
            end
        end
        ccov(ii,jj)=xmax;
    end
end
ccov=ccov+ccov';

% 计算近临函数矩阵lij
[x,ind]=sort(ccov);
m=size(ccov,1);
mij=zeros(m,m);
for ii=1:m
    for jj=1:m
        mij(ii,ind(jj,ii))=m-jj+1;
    end
end
lij=mij+mij'-2;
ee=6*eye(m);
lij=lij+ee;

% 构造初始聚类
[clu,rowid]=max(ccov);
[sclu,idc]=sort(clu);
sidc=idc;
cc=0;
while length(sidc)>0
    xx=rowid(sidc(1))*ones(1,length(sidc));
    xout=xx-sidc;
    if all(xout)
        cc=cc+1;
        catecell{cc}=sidc(1);
        sidc(1)=[];
    else
        cc=cc+1;
        catecell{cc}=[sidc(1) rowid(sidc(1))];
        for ii=2:length(sidc)
            if sidc(ii)==rowid(sidc(1));
                sidc(ii)=[];
                break
            end
        end
        sidc(1)=[];
    end
end

% 迭代聚类运算
dtime=0;
while length(catecell)>1
    % 计算类间最小近邻函数值rij
    dtime=dtime+1;
    rij=zeros(length(catecell),length(catecell));        
    for ii=1:length(catecell)-1
        ca=catecell{ii};
        for jj=ii+1:length(catecell)
            cb=catecell{jj};
            rij(ii,jj)=min(min(lij(ca,cb)));
        end
    end
    rij=rij+rij';
    rij=rij+diag(2*m*ones(size(rij,1),1));
    [ri,rowi]=min(rij);
    [sri,riind]=min(ri);
    bb=[];
    for ii=1:length(rowi)
        eli=catecell{ii};
        plij=lij(eli,eli);
        elong=length(eli);
        if elong>1
            cons=ones(elong,elong)-triu(ones(elong,elong));
            aimax=max(max(cons.*plij));
        else
            aimax=plij;
        end
        eli=catecell{rowi(ii)};
        plij=lij(eli,eli);
        elong=length(eli);
        if elong>1
            cons=ones(elong,elong)-triu(ones(elong,elong));
            akmax=max(max(cons.*plij));
        else
            akmax=plij;
        end
        if ri(ii)>aimax & ri(ii)>akmax
            bb(ii)=-((ri(ii)-aimax)+(ri(ii)-akmax));
        elseif ri(ii)<=aimax & ri(ii)>akmax
            bb(ii)=ri(ii)+aimax;
        elseif ri(ii)>aimax & ri(ii)<=akmax
            bb(ii)=ri(ii)+akmax;
        elseif ri(ii)<=aimax & ri(ii)<=akmax
            bb(ii)=ri(ii)+aimax+akmax;
        end
    end
    xjj(dtime)=sum(bb);
    if dtime==1,xjjmin=xjj(dtime);end
    if dtime>1 & xjj(dtime)<xjjmin
        xjjmin=xjj(dtime);
        clusterout=catecell;
    end
    if length(catecell)==2
        break
    end
     
    ieli=catecell{riind};
    plij=lij(ieli,ieli);
    pccov=ccov(ieli,ieli);
    elong=length(ieli);
    if elong>1
        cons=ones(elong,elong)-triu(ones(elong,elong));
        aimax=max(max(cons.*plij));
    else
        aimax=plij;
    end
    
    keli=catecell{rowi(riind)};
    plij=lij(keli,keli);
    elong=length(keli);
    if elong>1
        cons=ones(elong,elong)-triu(ones(elong,elong));
        akmax=max(max(cons.*plij));
    else
        akmax=plij;
    end
    
    if sri<=akmax | sri<=aimax
        sigmin=min(riind,rowi(riind));
        merg=catecell{sigmin};
        sigmax=max(riind,rowi(riind));
        merg=[merg catecell{sigmax}];
        catecell{sigmin}=merg;
        vcatecell=catecell;
        catecell=cell(1,length(vcatecell)-1);
        ccc=0;
        for ii=1:length(vcatecell);
            if ii-sigmax~=0
                ccc=ccc+1;
                catecell{ccc}=vcatecell{ii};
            end
        end
    else 
        clusterout=catecell;
        break
    end
end


% 计算聚类后各类的综合模板
for ii=1:length(clusterout) 
    cctv=clusterout{ii};
    clength=length(cctv);
    
    hmin=10000;wmin=10000;
    for kk=1:length(cctv);
        [hmod,wmod]=size(modstr(cctv(kk)).va);
        if hmod<hmin, hmin=hmod;end
        if wmod<wmin, wmin=wmod;end
    end
    
    ppccov=ccov(cctv,cctv);
    [vv,vidx]=max(ppccov);
    [v,idx]=max(vv);
    row=vidx(idx);col=idx;
    sump=sum(ppccov);
    if sump(vidx(idx))>=sump(idx)
        anum=cctv(vidx(idx));
        bnum=cctv(idx);
        a=modstr(cctv(vidx(idx))).va;
        b=modstr(cctv(idx)).va;
        [vx,ix]=sort(ppccov(row,:));
    else
        anum=cctv(idx);
        bnum=cctv(vidx(idx));
        a=modstr(cctv(idx)).va;
        b=modstr(cctv(vidx(idx))).va;
        [vx,ix]=sort(ppccov(col,:));
    end
    [ha,wa]=size(a);
    [hb,wb]=size(b);
    db=imdilate(b,sse);
    xmax=0;
    for kk=1:ha-hmin+1
        for ll=1:wa-wmin+1
            pa=a(kk:kk+hmin-1,ll:ll+wmin-1);
            spa=sum(sum(pa));
            for yy=1:hb-hmin+1
                for xx=1:wb-wmin+1
                    pdb=db(yy:yy+hmin-1,xx:xx+wmin-1);
                    out=and(pa,pdb);
                    percen=sum(sum(out))/spa;
                    if xmax<percen
                        xmax=percen;
                        outa=pa;
                        outb=b(yy:yy+hmin-1,xx:xx+wmin-1);
                        shifa=[kk,ll];
                        shifb=[yy,xx];
                        syout=out;
                    end
                end
            end
        end
    end
    leafstr=[];
    leafstr(1).va=outa;
    leafstr(1).sh=shifa;
    leafstr(1).su=sum(sum(outa));
    leafstr(1).nu=anum;
    
    leafstr(2).va=outb;
    leafstr(2).sh=shifb;
    leafstr(2).su=sum(sum(outb));
    leafstr(2).nu=bnum;
    cn=2;
    for i=length(ix)-1:-1:2
        sonmod=modstr(cctv(ix(i))).va;
        num=cctv(ix(i));
        dsonmod=imdilate(sonmod,sse);
        [ma,na]=size(dsonmod);
        xmax=0;
        sumb=sum(sum(syout));
        for kk=1:ma-hmin+1
            for ll=1:na-wmin+1
                outt=and(dsonmod(kk:kk+hmin-1,ll:ll+wmin-1),syout);
                percw=sum(sum(outt))/sumb;
                if xmax<percw;
                    xmax=percw;
                    omax=outt;
                    shif=[kk,ll];
                    pmod=sonmod(kk:kk+hmin-1,ll:ll+wmin-1);
                end
            end
        end
        cn=cn+1;
        syout=omax;
        leafstr(cn).va=pmod;
        leafstr(cn).sh=shif;
        leafstr(cn).su=sum(sum(pmod));
        leafstr(cn).nu=num;
    end
    sytree{1,1}=find(syout);
    sytree{1,2}=[hmin,wmin];
    sytree{2,1}=leafstr;
    syncell{ii}=sytree;
    sytree=[];
end

% matching
map=wwx([1400,2500,1000,2200]);
mapbw=edge(map,'canny');
dbw=imdilate(mapbw,sse);
obj=[];
[ih,iw]=size(map);
th=.85;flag=[];
for ii=1:length(syncell)
    sytree=syncell{ii};
    mothmod=sytree{1,1};
    hmin=sytree{1,2}(1);
    wmin=sytree{1,2}(2);
    leafstr=sytree{2,1};
    nums=length(leafstr);
    ssum=length(mothmod);
    for yy=1:ih-hmin+1
        for xx=1:iw-wmin+1
            pmap=dbw(yy:yy+hmin-1,xx:xx+wmin-1);
            percen=sum(pmap(mothmod))/ssum;
            if percen>=th
                flag(yy,xx)=1;
                pobj=[];
                for tt=1:length(leafstr)
                    percw=sum(sum(and(leafstr(tt).va,pmap)))/leafstr(tt).su;
                    if percw>th
                        pmapbw=mapbw(yy:yy+hmin-1,xx:xx+wmin-1);
                        dmod=imdilate(leafstr(tt).va,sse);
                        bpercw=sum(sum(and(pmapbw,dmod)))/sum(sum(pmapbw));
                        if bpercw>th
                            matcho=min(percw,bpercw);
                            if isempty(pobj)
                                row=yy-leafstr(tt).sh(1)+1;
                                col=xx-leafstr(tt).sh(2)+1;
                                pobj=[row,col,leafstr(tt).nu,matcho];
                            elseif matcho>pobj(4)
                                row=yy-leafstr(tt).sh(1)+1;
                                col=xx-leafstr(tt).sh(2)+1;
                                pobj=[row,col,leafstr(tt).nu,matcho];
                            end
                        end
                    end
                end
                obj=[obj;pobj];
            end
        end
    end
end

hmin=syncell{1}{1,2}(1);
wmin=syncell{1}{1,2}(2);
for i=2:length(syncell)
    if hmin>syncell{i}{1,2}(1)
        hmin=syncell{i}{1,2}(1);
    end
    if wmin>syncell{i}{1,2}(2)
        hmin=syncell{i}{1,2}(2);
    end
end
    
trueobj=[];
while length(obj)>0
    partobj=obj(1,:);
    flag=1;
    for i=2:length(obj)
        if abs(obj(i,1)-obj(1,1))<hmin & abs(obj(i,2)-obj(1,2))<wmin
            partobj=[partobj;obj(i,:)];
            flag=[flag i];
        end
    end
    [va,id]=max(partobj(:,4));
    trueobj=[trueobj;partobj(id,:)];
    obj(flag,:)=[];
end

toc