 Timestamp:
 Jan 18, 2020, 10:47:39 AM (4 years ago)
 Location:
 trunk/src
 Files:

 2 edited
Legend:
 Unmodified
 Added
 Removed

trunk/src/transform_field/ima2concentration.m
r1073 r1074 20 20 21 21 function [DataOut]=ima2concentration(DataIn,XmlData) 22 22 23 %% request input parameters 23 24 if isfield(DataIn,'Action') && isfield(DataIn.Action,'RUN') && isequal(DataIn.Action.RUN,0) … … 29 30 cpath=which('uvmat'); 30 31 addpath(fullfile(fileparts(cpath),'transform_field'))% define path for phys_polar.m 31 DataOut_1=[];32 32 33 34 %% for use in uvmat 35 % num_level=Data.ZIndex; 36 % if ~exist('Ref','var') 37 % huvmat=findobj(allchild(0),'tag','uvmat'); 38 % hhuvmat=guidata(huvmat); 39 % RootPath=get(hhuvmat.RootPath,'String'); 40 % 41 % %reference file 42 % RootPath=fullfile(RootPath,'LIF_REF'); 43 % file_ref=fullfile(RootPath,['lif_ref_' num2str(num_level) '.nc']); 44 % Ref=nc2struct(file_ref); 45 % end 46 47 %% Parameters 33 %% Transform images to polar coordinates with origin at the light source position 48 34 XmlData.TransformInput.PolarCentre=XmlData.LIFCalib.LightOrigin; %position of the laser origin [x, y] 49 50 % if isfield(XmlData.TransformInput,'PolarReferenceRadius')51 % def{2}=num2str(XmlData.TransformInput.PolarReferenceRadius);52 % end53 % if isfield(XmlData.TransformInput,'PolarReferenceAngle')54 % def{3}=num2str(XmlData.TransformInput.PolarReferenceAngle);55 56 %% concentration image57 35 DataIn.Action.RUN=1;% avoid input menu in phys_polar 58 36 DataOut=phys_polar(DataIn,XmlData); 59 % A=Ref.Aref;%default 60 % ind_good=find(Ref.Aref~=0); 61 % ind_bad=find(Ref.Aref==0); 62 % A(ind_good)=double(DataOut.A(ind_good))ImageOffset(1) 63 % %filtering and decimate 64 % Afilt=filter2(ones(nfilt,nfilt),A); 65 % Mask=filter2(ones(nfilt,nfilt),double(Ref.Aref~=0)); 66 % B=Afilt./Mask; 67 % A(ind_bad)=B(ind_bad); 68 % [npy,npx]=size(A); 69 % DataMask=DataOut; 70 % DataMask.A=2*ones(npy,npx);%mask=2 for good data 71 % 72 % DataMask.A(Ref.Aref==0)=1;%mask=0 for undefined data 73 % 74 % C=filter2(ones(nfilt,nfilt),Ref.Aref); 75 % D=C./Mask; 76 % Ref.Aref(ind_bad)=D(ind_bad); 77 % DataOut_1=[]; 78 % Coord_x=DataOut.Coord_x; 79 % Coord_y=DataOut.Coord_y; 80 % 81 % dX=(Coord_x(2)Coord_x(1))/(npx1); 82 % dY=(Coord_y(1)Coord_y(2))/(npy1);%mesh of new pixels 83 % [R,Y]=meshgrid(linspace(Coord_x(1),Coord_x(2),npx),linspace(Coord_y(1),Coord_y(2),npy)); 84 % r=Coord_x(1)+[0:npx1]*dX;%distance from laser 85 % %A(ind_good)=(A(ind_good)>=0).*A(ind_good); %replaces negative values by zeros 86 % A=A./Ref.Aref;% luminosity normalised by the reference (value at the edge of the box) 37 [npangle,npr]=size(DataOut.A);%size of the image in polar coordinates 38 dX=(DataOut.Coord_x(2)DataOut.Coord_x(1))/(npr1);% radial step 87 39 88 [npangle,npr]=size(DataOut.A); 89 dX=(DataOut.Coord_x(2)DataOut.Coord_x(1))/(npr1); 40 %% introduce the reference line where the laser enters the fluid region 90 41 r_edge=XmlData.LIFCalib.RefLineRadius'*ones(1,npr);% radial position of the reference line extended as a matrix (npx,npy) 91 A_ref=XmlData.LIFCalib.RefLineLum'*ones(1,npr);% luminosity on the reference line at the edge of the box,extended as a matrix (npx,npy)42 A_ref=XmlData.LIFCalib.RefLineLum'*ones(1,npr);% luminosity on the reference line extended as a matrix (npx,npy) 92 43 R=ones(npangle,1)*linspace(DataOut.Coord_x(1), DataOut.Coord_x(2),npr);%radial coordinate extended as a matrix (npx,npy) 93 %Edge_ind=find((abs(Rr_edge)/dX)<=1 & DataMask.A~=0);%indies of positions close to r_edge, values greater than 1 are not expected94 %yedge=min(min(Y(Edge_ind)));95 % jmax=round((yedgeCoord_y(1))/dY+1);96 % DataMask.A(jmax:end,:)=0;97 %98 % A(isnan(A)isinf(A))=0;99 44 100 % radius along the reference line 101 %Theta=(linspace(Coord_y(1),Coord_y(2),npy)*pi/180)'*ones(1,npx);%theta in radians 102 103 gamma_coeff=XmlData.LIFCalib.DecayRate; 104 105 45 %gamma_coeff=XmlData.LIFCalib.DecayRate; 106 46 DataOut.A(R<r_edge)=0; 107 DataOut.A=double(DataOut.A)./A_ref; 108 I=(r_edgedX* gamma_coeff.*cumsum(R.*DataOut.A,2))./R;% expected laser intensity along the line109 DataOut.A=DataOut.A./I;%concentration 47 DataOut.A=double(DataOut.A)./A_ref;% renormalize the luminosity with the reference luminosity at the same azimuth on the reference line 48 I=(r_edgedX*XmlData.LIFCalib.DecayRate.*cumsum(R.*DataOut.A,2))./R;% expected laser intensity along the line 49 DataOut.A=DataOut.A./I;%concentration normalized by the uniform concentration assumed in the ref image used for calibration 110 50 DataOut.A(I<=0)=0;% eliminate values obtained with I<=0 111 51 112 RangeX=DataIn.Coord_xXmlData.LIFCalib.LightOrigin(1); 113 RangeY=DataIn.Coord_yXmlData.LIFCalib.LightOrigin(2); 114 % 115 DataOut=polar2phys(DataOut,RangeX,RangeY); 116 DataOut.A=uint16(DataOut.A); 117 DataOut.Coord_x=DataOut.Coord_x+XmlData.LIFCalib.LightOrigin(1); 52 DataOut=polar2phys(DataOut);% back to phys cartesian coordinates with origin at the light source 53 DataOut.A=uint16(1000*DataOut.A);% concentration multiplied by 1000 to get an image 54 DataOut.Coord_x=DataOut.Coord_x+XmlData.LIFCalib.LightOrigin(1);%shift to original cartesian coordinates 118 55 DataOut.Coord_y=DataOut.Coord_y+XmlData.LIFCalib.LightOrigin(2); 119 56 120 57 58 function DataOut=polar2phys(DataIn) 59 %%%%%%%%%%%%%%%%%%%% 60 DataOut=DataIn; %default 61 [npy,npx]=size(DataIn.A); 62 dx=(DataIn.Coord_x(2)DataIn.Coord_x(1))/(npx1); %mesh along radius 63 dy=(DataIn.Coord_y(2)DataIn.Coord_y(1))/(npy1);%mesh along azimuth 121 64 122 function DataOut=polar2phys(DataIn,RangeX,RangeY) 123 %%%%%%%%%%%%%%%%%%%% 124 DataOut=DataIn; %fdefault 125 [npy,npx]=size(DataIn.A); 126 dx=(DataIn.Coord_x(2)DataIn.Coord_x(1))/(npx1); 127 dy=(DataIn.Coord_y(2)DataIn.Coord_y(1))/(npy1);%mesh 65 %% create cartesian coordinates in the domain defined by the four image corners 128 66 rcorner=[DataIn.Coord_x(1) DataIn.Coord_x(2) DataIn.Coord_x(1) DataIn.Coord_x(2)];% radius of the corners 129 67 ycorner=[DataIn.Coord_y(2) DataIn.Coord_y(2) DataIn.Coord_y(1) DataIn.Coord_y(1)];% azimuth of the corners 130 68 thetacorner=pi*ycorner/180;% azimuth in radians 131 69 [Xcorner,Ycorner] = pol2cart(thetacorner,rcorner);% cartesian coordinates of the corners (with respect to lser source) 132 if ~exist('RangeX','var')133 70 RangeX(1)=min(Xcorner); 134 71 RangeX(2)=max(Xcorner); 135 end136 if ~exist('RangeY','var')137 72 RangeY(2)=min(Ycorner); 138 73 RangeY(1)=max(Ycorner); 139 end140 %Rangx=[100 100];%bounds of the initial box141 %Rangy=[75 150];142 % Rangy(1)=min(Ycorner);143 % Rangy(2)=max(Ycorner);144 74 x=linspace(RangeX(1),RangeX(2),npx);%coordinates of the new pixels 145 75 y=linspace(RangeY(2),RangeY(1),npy); 146 [X,Y]=meshgrid(x,y);%grid for new pixels in cartesian coordi antes76 [X,Y]=meshgrid(x,y);%grid for new pixels in cartesian coordinates 147 77 78 %% image indices corresponding to the cartesian grid 148 79 [Theta,R] = cart2pol(X,Y);%corresponding polar coordiantes 149 Theta=Theta*180/pi; 150 %Theta=1+round((ThetaDataIn.Coord_y(1))/dy); %index along y (dy negative) 151 Theta=1round((ThetaDataIn.Coord_y(2))/dy); %index along y (dy negative) 152 R=1+round((RDataIn.Coord_x(1))/dx); %index along x 80 Theta=180*Theta/pi;%angles in degrees 81 Theta=1round((ThetaDataIn.Coord_y(2))/dy); %angular index along y (dy negative) 82 R=1+round((RDataIn.Coord_x(1))/dx); %angular index along x 153 83 R=reshape(R,1,npx*npy);%indices reorganized in 'line' 154 84 Theta=reshape(Theta,1,npx*npy); … … 162 92 vec_B(ind_out)=zeros(size(ind_out)); 163 93 DataOut.A=flipdim(reshape(vec_B,npy,npx),1);%new image in real coordinates 164 165 %Rangx=Rangxradius_ref;166 94 DataOut.Coord_x=RangeX; 167 95 DataOut.Coord_y=RangeY; 
trunk/src/uvmat.m
r1072 r1074 1439 1439 select_line=zeros(1,numel(ListObj)); 1440 1440 select_mask=zeros(1,numel(ListObj)); 1441 for iobj=1:numel(ListObj) ;1442 if isfield(ListObj{iobj},'Type') && strcmp(ListObj{iobj}.Type,'line') 1441 for iobj=1:numel(ListObj) 1442 if isfield(ListObj{iobj},'Type') && strcmp(ListObj{iobj}.Type,'line')&& ~strcmp(ListObj{iobj}.Name,'MaskPolygon') 1443 1443 select_line(iobj)=1;% select the lines among the projection objects 1444 1444 end
Note: See TracChangeset
for help on using the changeset viewer.