Ignore:
Timestamp:
Jun 11, 2016, 9:32:29 PM (8 years ago)
Author:
sommeria
Message:

various updates

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/series/sub_background.m

    r924 r950  
    340340        if Param.ActionInput.CheckLevelTransform
    341341            C=levels(Acor);
    342             imwrite(C,newname,'BitDepth',8); % save the new image
     342            imwrite(C,newname,'BitDepth',16); % save the new image
    343343        else
    344344            if isequal(FileInfo{1}.BitDepth,16)
     
    393393                if Param.ActionInput.CheckLevelTransform
    394394                    C=levels(Acor);
    395                     imwrite(C,newname,'BitDepth',8); % save the new image
     395                    imwrite(C,newname,'BitDepth',16); % save the new image
    396396                else
    397397                    if isequal(FileInfo{1}.BitDepth,16)
     
    421421        if Param.ActionInput.CheckLevelTransform
    422422            C=levels(Acor);
    423             imwrite(C,newname,'BitDepth',8); % save the new image
     423            imwrite(C,newname,'BitDepth',16); % save the new image
    424424        else
    425425            if isequal(FileInfo{1}.BitDepth,16)
     
    435435end
    436436
    437 
    438437function C=levels(A)
    439 %whos A;
    440 B=double(A(:,:,1));
    441 windowsize=round(min(size(B,1),size(B,2))/20);
    442 windowsize=floor(windowsize/2)*2+1;
    443 ix=1/2-windowsize/2:-1/2+windowsize/2;%
    444 %del=np/3;
    445 %fct=exp(-(ix/del).^2);
    446 fct2=cos(ix/(windowsize-1)/2*pi/2);
    447 %Mfiltre=(ones(5,5)/5^2);
    448 %Mfiltre=fct2';
    449 Mfiltre=fct2'*fct2;
    450 Mfiltre=Mfiltre/(sum(sum(Mfiltre)));
    451 
    452 C=filter2(Mfiltre,B);
    453 C(:,1:windowsize)=C(:,windowsize)*ones(1,windowsize);
    454 C(:,end-windowsize+1:end)=C(:,end-windowsize+1)*ones(1,windowsize);
    455 C(1:windowsize,:)=ones(windowsize,1)*C(windowsize,:);
    456 C(end-windowsize+1:end,:)=ones(windowsize,1)*C(end-windowsize,:);
    457 C=tanh(B./(2*C));
    458 [n,c]=hist(reshape(C,1,[]),100);
    459 % figure;plot(c,n);
    460 
    461 [m,i]=max(n);
    462 c_max=c(i);
    463 [dummy,index]=sort(abs(c-c(i)));
    464 n=n(index);
    465 c=c(index);
    466 i_select = find(cumsum(n)<0.95*sum(n));
    467 if isempty(i_select)
    468     i_select = 1:length(c);
    469 end
    470 c_select=c(i_select);
    471 n_select=n(i_select);
    472 cmin=min(c_select);
    473 cmax=max(c_select);
    474 C=(C-cmin)/(cmax-cmin)*256;
    475 C=uint8(C);
     438
     439nblock_y=100;%2*Param.TransformInput.BlockSize;
     440nblock_x=100;%2*Param.TransformInput.BlockSize;
     441[npy,npx]=size(A);
     442[X,Y]=meshgrid(1:npx,1:npy);
     443
     444%Backg=zeros(size(A));
     445%Aflagmin=sparse(imregionalmin(A));%Amin=1 for local image minima
     446%Amin=A.*Aflagmin;%values of A at local minima
     447% local background: find all the local minima in image subblocks
     448fctblock= inline('median(x(:))');
     449Backg=blkproc(A,[nblock_y nblock_x],fctblock);% take the median in  blocks
     450fctblock= inline('mean(x(:))');
     451B=imresize(Backg,size(A),'bilinear');% interpolate to the initial size image
     452A=(A-B);%substract background
     453AMean=blkproc(A,[nblock_y nblock_x],fctblock);% take the mean in  blocks
     454fctblock= inline('var(x(:))');
     455AVar=blkproc(A,[nblock_y nblock_x],fctblock);% take the mean in  blocks
     456Avalue=AVar./AMean% typical value of particle luminosity
     457Avalue=imresize(Avalue,size(A),'bilinear');% interpolate to the initial size image
     458C=uint16(1000*tanh(A./(2*Avalue)));
     459%Bmin=blkproc(Aflagmin,[nblock_y nblock_x],sumblock);% find the number of minima in blocks
     460%Backg=Backg./Bmin; % find the average of minima in blocks
     461% function C=levels(A)
     462% %whos A;
     463% B=double(A(:,:,1));
     464% windowsize=round(min(size(B,1),size(B,2))/20);
     465% windowsize=floor(windowsize/2)*2+1;
     466% ix=1/2-windowsize/2:-1/2+windowsize/2;%
     467% %del=np/3;
     468% %fct=exp(-(ix/del).^2);
     469% fct2=cos(ix/(windowsize-1)/2*pi/2);
     470% %Mfiltre=(ones(5,5)/5^2);
     471% %Mfiltre=fct2';
     472% Mfiltre=fct2'*fct2;
     473% Mfiltre=Mfiltre/(sum(sum(Mfiltre)));
     474%
     475% C=filter2(Mfiltre,B);
     476% C(:,1:windowsize)=C(:,windowsize)*ones(1,windowsize);
     477% C(:,end-windowsize+1:end)=C(:,end-windowsize+1)*ones(1,windowsize);
     478% C(1:windowsize,:)=ones(windowsize,1)*C(windowsize,:);
     479% C(end-windowsize+1:end,:)=ones(windowsize,1)*C(end-windowsize,:);
     480% C=tanh(B./(2*C));
     481% [n,c]=hist(reshape(C,1,[]),100);
     482% % figure;plot(c,n);
     483%
     484% [m,i]=max(n);
     485% c_max=c(i);
     486% [dummy,index]=sort(abs(c-c(i)));
     487% n=n(index);
     488% c=c(index);
     489% i_select = find(cumsum(n)<0.95*sum(n));
     490% if isempty(i_select)
     491%     i_select = 1:length(c);
     492% end
     493% c_select=c(i_select);
     494% n_select=n(i_select);
     495% cmin=min(c_select);
     496% cmax=max(c_select);
     497% C=(C-cmin)/(cmax-cmin)*256;
     498% C=uint8(C);
Note: See TracChangeset for help on using the changeset viewer.