source: trunk/src/uvmat_doc/FUNCTIONS_DOC/phys_polar.html @ 37

Last change on this file since 37 was 37, checked in by sommeria, 14 years ago

create_grid.fig ,
uvmat_doc and all the included files added
File size: 19.7 KB
Line 
1<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN"
2                "http://www.w3.org/TR/REC-html40/loose.dtd">
3<html>
4<head>
5  <title>Description of phys_polar</title>
6  <meta name="keywords" content="phys_polar">
7  <meta name="description" content="transform image coordinates (px) to physical coordinates">
8  <meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1">
9  <meta name="generator" content="m2html &copy; 2003 Guillaume Flandin">
10  <meta name="robots" content="index, follow">
11  <link type="text/css" rel="stylesheet" href="../m2html.css">
12</head>
13<body>
14<a name="_top"></a>
15<div><a href="../index.html">Home</a> &gt;  <a href="index.html">.</a> &gt; phys_polar.m</div>
16
17<!--<table width="100%"><tr><td align="left"><a href="../index.html"><img alt="<" border="0" src="../left.png">&nbsp;Master index</a></td>
18<td align="right"><a href="index.html">Index for .&nbsp;<img alt=">" border="0" src="../right.png"></a></td></tr></table>-->
19
20<h1>phys_polar
21</h1>
22
23<h2><a name="_name"></a>PURPOSE <a href="#_top"><img alt="^" border="0" src="../up.png"></a></h2>
24<div class="box"><strong>transform image coordinates (px) to physical coordinates</strong></div>
25
26<h2><a name="_synopsis"></a>SYNOPSIS <a href="#_top"><img alt="^" border="0" src="../up.png"></a></h2>
27<div class="box"><strong>function [DataOut,DataOut_1]=phys_polar(varargin) </strong></div>
28
29<h2><a name="_description"></a>DESCRIPTION <a href="#_top"><img alt="^" border="0" src="../up.png"></a></h2>
30<div class="fragment"><pre class="comment">transform image coordinates (px) to physical coordinates
31 then transform to polar coordinates:
32[DataOut,DataOut_1]=phys_polar(varargin)
33
34 OUTPUT:
35 DataOut: structure of modified data field: .X=radius, .Y=azimuth angle, .U, .V are radial and azimuthal velocity components
36 DataOut_1:  second data field (if two fields are in input)
37
38INPUT:
39 Data:  structure of input data (like UvData)
40 CalibData= structure containing the field .GeometryCalib with calibration parameters
41 Data_1:  second input field (not mandatory)
42 CalibData_1= calibration parameters for the second field</pre></div>
43
44<!-- crossreference -->
45<h2><a name="_cross"></a>CROSS-REFERENCE INFORMATION <a href="#_top"><img alt="^" border="0" src="../up.png"></a></h2>
46This function calls:
47<ul style="list-style-image:url(../matlabicon.gif)">
48<li><a href="phys_XYZ.html" class="code" title="function [Xphys,Yphys,Zphys]=phys_XYZ(Calib,X,Y,Z)">phys_XYZ</a>        'phys_XYZ':transforms image (px) to real world (phys) coordinates using geometric calibration parameters</li><li><a href="px_XYZ.html" class="code" title="function [X,Y]=px_XYZ(Calib,Xphys,Yphys,Zphys)">px_XYZ</a>   'px_XYZ': transform phys coordinates to image coordinates (px)</li></ul>
49This function is called by:
50<ul style="list-style-image:url(../matlabicon.gif)">
51</ul>
52<!-- crossreference -->
53
54<h2><a name="_subfunctions"></a>SUBFUNCTIONS <a href="#_top"><img alt="^" border="0" src="../up.png"></a></h2>
55<ul style="list-style-image:url(../matlabicon.gif)">
56<li><a href="#_sub1" class="code">function DataOut=phys_1(Data,Calib,origin_xy,radius_offset,angle_offset,angle_scale)</a></li><li><a href="#_sub2" class="code">function [A_out,Rangx,Rangy]=phys_Ima(A,CalibIn,ZIndex,origin_xy,radius_offset,angle_offset,angle_scale)</a></li></ul>
57<h2><a name="_source"></a>SOURCE CODE <a href="#_top"><img alt="^" border="0" src="../up.png"></a></h2>
58<div class="fragment"><pre>0001 <span class="comment">%transform image coordinates (px) to physical coordinates</span>
590002 <span class="comment">% then transform to polar coordinates:</span>
600003 <span class="comment">%[DataOut,DataOut_1]=phys_polar(varargin)</span>
610004 <span class="comment">%</span>
620005 <span class="comment">% OUTPUT:</span>
630006 <span class="comment">% DataOut: structure of modified data field: .X=radius, .Y=azimuth angle, .U, .V are radial and azimuthal velocity components</span>
640007 <span class="comment">% DataOut_1:  second data field (if two fields are in input)</span>
650008 <span class="comment">%</span>
660009 <span class="comment">%INPUT:</span>
670010 <span class="comment">% Data:  structure of input data (like UvData)</span>
680011 <span class="comment">% CalibData= structure containing the field .GeometryCalib with calibration parameters</span>
690012 <span class="comment">% Data_1:  second input field (not mandatory)</span>
700013 <span class="comment">% CalibData_1= calibration parameters for the second field</span>
710014
720015 <a name="_sub0" href="#_subfunctions" class="code">function [DataOut,DataOut_1]=phys_polar(varargin)</a>
730016 Calib{1}=[];
740017 <span class="keyword">if</span> nargin==2||nargin==4
750018     Data=varargin{1};
760019     DataOut=Data;<span class="comment">%default</span>
770020     DataOut_1=[];<span class="comment">%default</span>
780021     CalibData=varargin{2};
790022     <span class="keyword">if</span> isfield(CalibData,<span class="string">'GeometryCalib'</span>)
800023         Calib{1}=CalibData.GeometryCalib;
810024     <span class="keyword">end</span>
820025     Calib{2}=Calib{1};
830026 <span class="keyword">else</span>
840027     DataOut.Txt=<span class="string">'wrong input: need two or four structures'</span>;
850028 <span class="keyword">end</span>
860029 test_1=0;
870030 <span class="keyword">if</span> nargin==4
880031     test_1=1;
890032     Data_1=varargin{3};
900033     DataOut_1=Data_1;<span class="comment">%default</span>
910034     CalibData_1=varargin{4};
920035     <span class="keyword">if</span> isfield(CalibData_1,<span class="string">'GeometryCalib'</span>)
930036         Calib{2}=CalibData_1.GeometryCalib;
940037     <span class="keyword">end</span>
950038 <span class="keyword">end</span>
960039
970040 <span class="comment">%parameters for polar coordinates (taken from the calibration data of the first field)</span>
980041 <span class="comment">%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span>
990042 origin_xy=[0 0];<span class="comment">%center for the polar coordinates in the original x,y coordinates</span>
1000043 <span class="keyword">if</span> isfield(Calib{1},<span class="string">'PolarCentre'</span>) &amp;&amp; isnumeric(Calib{1}.PolarCentre)
1010044     <span class="keyword">if</span> isequal(length(Calib{1}.PolarCentre),2);
1020045         origin_xy= Calib{1}.PolarCentre;
1030046     <span class="keyword">end</span>
1040047 <span class="keyword">end</span>
1050048 radius_offset=0;<span class="comment">%reference radius used to offset the radial coordinate r</span>
1060049 angle_offset=0; <span class="comment">%reference angle used as new origin of the polar angle (= axis Ox by default)</span>
1070050 <span class="keyword">if</span> isfield(Calib{1},<span class="string">'PolarReferenceRadius'</span>) &amp;&amp; isnumeric(Calib{1}.PolarReferenceRadius)
1080051     radius_offset=Calib{1}.PolarReferenceRadius;
1090052 <span class="keyword">end</span>
1100053 <span class="keyword">if</span> radius_offset &gt; 0
1110054     angle_scale=radius_offset; <span class="comment">%the azimuth is rescale in terms of the length along the reference radius</span>
1120055 <span class="keyword">else</span>
1130056     angle_scale=180/pi; <span class="comment">%polar angle in degrees</span>
1140057 <span class="keyword">end</span>
1150058 <span class="keyword">if</span> isfield(Calib{1},<span class="string">'PolarReferenceAngle'</span>) &amp;&amp; isnumeric(Calib{1}.PolarReferenceAngle)
1160059     angle_offset=Calib{1}.PolarReferenceAngle; <span class="comment">%offset angle (in unit of the final angle, degrees or arc length along the reference radius))</span>
1170060 <span class="keyword">end</span>
1180061 <span class="comment">% new x coordinate = radius-radius_offset;</span>
1190062 <span class="comment">% new y coordinate = theta*angle_scale-angle_offset</span>
1200063
1210064 <span class="comment">%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span>
1220065
1230066 iscalar=0;
1240067 <span class="keyword">if</span>  ~isempty(Calib{1})
1250068     DataOut=<a href="#_sub1" class="code" title="subfunction DataOut=phys_1(Data,Calib,origin_xy,radius_offset,angle_offset,angle_scale)">phys_1</a>(Data,Calib{1},origin_xy,radius_offset,angle_offset,angle_scale);
1260069     <span class="comment">%case of images or scalar</span>
1270070     <span class="keyword">if</span> isfield(Data,<span class="string">'A'</span>)&amp;isfield(Data,<span class="string">'AX'</span>)&amp;~isempty(Data.AX) &amp; isfield(Data,<span class="string">'AY'</span>)&amp;<span class="keyword">...</span>
1280071                                            ~isempty(Data.AY)&amp;length(Data.A)&gt;1
1290072         iscalar=1;
1300073         A{1}=Data.A;
1310074     <span class="keyword">end</span>
1320075     <span class="comment">%transform of X,Y coordinates for vector fields</span>
1330076     <span class="keyword">if</span> isfield(Data,<span class="string">'ZIndex'</span>)&amp;~isempty(Data.ZIndex)
1340077         ZIndex=Data.ZIndex;
1350078     <span class="keyword">else</span>
1360079         ZIndex=0;
1370080     <span class="keyword">end</span>
1380081 <span class="keyword">end</span>
1390082
1400083 <span class="keyword">if</span> test_1
1410084     DataOut_1=<a href="#_sub1" class="code" title="subfunction DataOut=phys_1(Data,Calib,origin_xy,radius_offset,angle_offset,angle_scale)">phys_1</a>(Data_1,Calib{2},origin_xy,radius_offset,angle_offset,angle_scale);
1420085     <span class="keyword">if</span> isfield(Data_1,<span class="string">'A'</span>)&amp;isfield(Data_1,<span class="string">'AX'</span>)&amp;~isempty(Data_1.AX) &amp; isfield(Data_1,<span class="string">'AY'</span>)&amp;<span class="keyword">...</span>
1430086                                        ~isempty(Data_1.AY)&amp;length(Data_1.A)&gt;1
1440087           iscalar=iscalar+1;
1450088           Calib{iscalar}=Calib{2};
1460089           A{iscalar}=Data_1.A;
1470090           <span class="keyword">if</span> isfield(Data_1,<span class="string">'ZIndex'</span>)&amp;~isequal(Data_1.ZIndex,ZIndex)
1480091               DataOut.Txt=<span class="string">'inconsistent plane indexes in the two input fields'</span>;
1490092           <span class="keyword">end</span>
1500093           <span class="keyword">if</span> iscalar==1<span class="comment">% case for which only the second field is a scalar</span>
1510094                [A,AX,AY]=<a href="#_sub2" class="code" title="subfunction [A_out,Rangx,Rangy]=phys_Ima(A,CalibIn,ZIndex,origin_xy,radius_offset,angle_offset,angle_scale)">phys_Ima</a>(A,Calib,ZIndex,origin_xy,radius_offset,angle_offset,angle_scale);
1520095                DataOut_1.A=A{1};
1530096                DataOut_1.AX=AX;
1540097                DataOut_1.AY=AY;
1550098                <span class="keyword">return</span>
1560099           <span class="keyword">end</span>
1570100     <span class="keyword">end</span>
1580101 <span class="keyword">end</span>
1590102 <span class="keyword">if</span> iscalar~=0
1600103     [A,AX,AY]=<a href="#_sub2" class="code" title="subfunction [A_out,Rangx,Rangy]=phys_Ima(A,CalibIn,ZIndex,origin_xy,radius_offset,angle_offset,angle_scale)">phys_Ima</a>(A,Calib,ZIndex,origin_xy,radius_offset,angle_offset,angle_scale);<span class="comment">%</span>
1610104     DataOut.A=A{1};
1620105     DataOut.AX=AX;
1630106     DataOut.AY=AY;
1640107     <span class="keyword">if</span> iscalar==2
1650108         DataOut_1.A=A{2};
1660109         DataOut_1.AX=AX;
1670110         DataOut_1.AY=AY;
1680111     <span class="keyword">end</span>
1690112 <span class="keyword">end</span>
1700113
1710114 <span class="comment">%------------------------------------------------</span>
1720115 <a name="_sub1" href="#_subfunctions" class="code">function DataOut=phys_1(Data,Calib,origin_xy,radius_offset,angle_offset,angle_scale)</a>
1730116
1740117 DataOut=Data;
1750118 DataOut.CoordType=<span class="string">'phys'</span>; <span class="comment">%put flag for physical coordinates</span>
1760119 <span class="keyword">if</span> isfield(Calib,<span class="string">'CoordUnit'</span>)
1770120     DataOut.CoordUnit=Calib.CoordUnit;
1780121 <span class="keyword">else</span>
1790122     DataOut.CoordUnit=<span class="string">'cm'</span>; <span class="comment">%default</span>
1800123 <span class="keyword">end</span>
1810124 DataOut.TimeUnit=<span class="string">'s'</span>;
1820125 <span class="comment">%perform a geometry transform if Calib contains a field .GeometryCalib</span>
1830126 <span class="keyword">if</span> isfield(Data,<span class="string">'CoordType'</span>) &amp;&amp; isequal(Data.CoordType,<span class="string">'px'</span>) &amp;&amp; ~isempty(Calib)
1840127     <span class="keyword">if</span> isfield(Data,<span class="string">'CoordUnit'</span>)
1850128          DataOut=rmfield(DataOut,<span class="string">'CoordUnit'</span>);
1860129     <span class="keyword">end</span>
1870130     <span class="comment">%transform of X,Y coordinates for vector fields</span>
1880131     <span class="keyword">if</span> isfield(Data,<span class="string">'ZIndex'</span>)&amp;~isempty(Data.ZIndex)
1890132         Z=Data.ZIndex;
1900133     <span class="keyword">else</span>
1910134         Z=0;
1920135     <span class="keyword">end</span>
1930136     <span class="keyword">if</span> isfield(Data,<span class="string">'X'</span>) &amp;isfield(Data,<span class="string">'Y'</span>)&amp;~isempty(Data.X) &amp; ~isempty(Data.Y)
1940137         [DataOut.X,DataOut.Y,DataOut.Z]=<a href="phys_XYZ.html" class="code" title="function [Xphys,Yphys,Zphys]=phys_XYZ(Calib,X,Y,Z)">phys_XYZ</a>(Calib,Data.X,Data.Y,Z); <span class="comment">%transform from pixels to physical</span>
1950138         DataOut.X=DataOut.X-origin_xy(1);<span class="comment">%origin of coordinates at the tank center</span>
1960139         DataOut.Y=DataOut.Y-origin_xy(2);<span class="comment">%origin of coordinates at the tank center</span>
1970140         [theta,DataOut.X] = cart2pol(DataOut.X,DataOut.Y);<span class="comment">%theta  and X are the polar coordinates angle and radius</span>
1980141           <span class="comment">%shift and renormalize the polar coordinates</span>
1990142         DataOut.X=DataOut.X-radius_offset;<span class="comment">%</span>
2000143         DataOut.Y=theta*angle_scale-angle_offset;<span class="comment">% normalized angle: distance along reference radius</span>
2010144         <span class="comment">%transform velocity field if exists</span>
2020145         <span class="keyword">if</span> isfield(Data,<span class="string">'U'</span>)&amp;isfield(Data,<span class="string">'V'</span>)&amp;~isempty(Data.U) &amp; ~isempty(Data.V)&amp; isfield(Data,<span class="string">'dt'</span>)
2030146             <span class="keyword">if</span> ~isempty(Data.dt)
2040147             [XOut_1,YOut_1]=<a href="phys_XYZ.html" class="code" title="function [Xphys,Yphys,Zphys]=phys_XYZ(Calib,X,Y,Z)">phys_XYZ</a>(Calib,Data.X-Data.U/2,Data.Y-Data.V/2,Z);
2050148             [XOut_2,YOut_2]=<a href="phys_XYZ.html" class="code" title="function [Xphys,Yphys,Zphys]=phys_XYZ(Calib,X,Y,Z)">phys_XYZ</a>(Calib,Data.X+Data.U/2,Data.Y+Data.V/2,Z);
2060149             UX=(XOut_2-XOut_1)/Data.dt;
2070150             VY=(YOut_2-YOut_1)/Data.dt;     
2080151             <span class="comment">%transform u,v into polar coordiantes</span>
2090152             DataOut.U=UX.*cos(theta)+VY.*sin(theta);<span class="comment">%radial velocity</span>
2100153             DataOut.V=(-UX.*sin(theta)+VY.*cos(theta));<span class="comment">%./(DataOut.X)%+radius_ref);%angular velocity calculated</span>
2110154             <span class="comment">%shift and renormalize the angular velocity</span>
2120155             <span class="keyword">end</span>
2130156         <span class="keyword">end</span>
2140157     <span class="keyword">end</span>
2150158 <span class="keyword">end</span>
2160159
2170160 
2180161 <span class="comment">%%%%%%%%%%%%%%%%%%%%</span>
2190162 <a name="_sub2" href="#_subfunctions" class="code">function [A_out,Rangx,Rangy]=phys_Ima(A,CalibIn,ZIndex,origin_xy,radius_offset,angle_offset,angle_scale)</a>
2200163 xcorner=[];
2210164 ycorner=[];
2220165 npx=[];
2230166 npy=[];
2240167
2250168 <span class="keyword">for</span> icell=1:length(A)
2260169     siz=size(A{icell});
2270170     npx=[npx siz(2)];
2280171     npy=[npy siz(1)];
2290172     zphys=0; <span class="comment">%default</span>
2300173     <span class="keyword">if</span> isfield(CalibIn{icell},<span class="string">'SliceCoord'</span>) <span class="comment">%.Z= index of plane</span>
2310174        SliceCoord=CalibIn{icell}.SliceCoord(ZIndex,:);
2320175        zphys=SliceCoord(3); <span class="comment">%to generalize for non-parallel planes</span>
2330176     <span class="keyword">end</span>
2340177     xima=[0.5 siz(2)-0.5 0.5 siz(2)-0.5];<span class="comment">%image coordiantes of corners</span>
2350178     yima=[0.5 0.5 siz(1)-0.5 siz(1)-0.5];
2360179     [xcorner_new,ycorner_new]=<a href="phys_XYZ.html" class="code" title="function [Xphys,Yphys,Zphys]=phys_XYZ(Calib,X,Y,Z)">phys_XYZ</a>(CalibIn{icell},xima,yima,ZIndex);<span class="comment">%corresponding physical coordinates</span>
2370180     <span class="comment">%transform the corner coordinates into polar ones</span>
2380181     xcorner_new=xcorner_new-origin_xy(1);<span class="comment">%shift to the origin of the polar coordinates</span>
2390182     ycorner_new=ycorner_new-origin_xy(2);<span class="comment">%shift to the origin of the polar coordinates</span>
2400183     [theta,xcorner_new] = cart2pol(xcorner_new,ycorner_new);<span class="comment">%theta  and X are the polar coordinates angle and radius</span>
2410184     <span class="keyword">if</span> (max(theta)-min(theta))&gt;pi   <span class="comment">%if the polar origin is inside the image</span>
2420185         xcorner_new=[0 max(xcorner_new)];
2430186         theta=[-pi pi];
2440187     <span class="keyword">end</span>
2450188           <span class="comment">%shift and renormalize the polar coordinates</span>
2460189     xcorner_new=xcorner_new-radius_offset;<span class="comment">%</span>
2470190     ycorner_new=theta*angle_scale-angle_offset;<span class="comment">% normalized angle: distance along reference radius</span>
2480191     xcorner=[xcorner xcorner_new];
2490192     ycorner=[ycorner ycorner_new];
2500193 <span class="keyword">end</span>
2510194 Rangx(1)=min(xcorner);
2520195 Rangx(2)=max(xcorner);
2530196 Rangy(2)=min(ycorner);
2540197 Rangy(1)=max(ycorner);
2550198 <span class="comment">% test_multi=(max(npx)~=min(npx)) | (max(npy)~=min(npy));</span>
2560199 npx=max(npx);
2570200 npy=max(npy);
2580201 x=linspace(Rangx(1),Rangx(2),npx);
2590202 y=linspace(Rangy(1),Rangy(2),npy);
2600203 [X,Y]=meshgrid(x,y);<span class="comment">%grid in physical coordinates</span>
2610204 <span class="comment">%transform X, Y in cartesian</span>
2620205 X=X+radius_offset;<span class="comment">%</span>
2630206 Y=(Y+angle_offset)/angle_scale;<span class="comment">% normalized angle: distance along reference radius</span>
2640207 [X,Y] = pol2cart(Y,X);
2650208 X=X+origin_xy(1);<span class="comment">%shift to the origin of the polar coordinates</span>
2660209 Y=Y+origin_xy(2);<span class="comment">%shift to the origin of the polar coordinates</span>
2670210 <span class="keyword">for</span> icell=1:length(A)
2680211     [XIMA,YIMA]=<a href="px_XYZ.html" class="code" title="function [X,Y]=px_XYZ(Calib,Xphys,Yphys,Zphys)">px_XYZ</a>(CalibIn{icell},X,Y,zphys);<span class="comment">%corresponding image indices for each point in the real space grid</span>
2690212     XIMA=reshape(round(XIMA),1,npx*npy);<span class="comment">%indices reorganized in 'line'</span>
2700213     YIMA=reshape(round(YIMA),1,npx*npy);
2710214     flagin=XIMA&gt;=1 &amp; XIMA&lt;=npx &amp; YIMA &gt;=1 &amp; YIMA&lt;=npy;<span class="comment">%flagin=1 inside the original image</span>
2720215     vec_A=reshape(A{icell}(:,:,1),1,npx*npy);<span class="comment">%put the original image in line</span>
2730216     ind_in=find(flagin);
2740217     ind_out=find(~flagin);
2750218     ICOMB=((XIMA-1)*npy+(npy+1-YIMA));
2760219     ICOMB=ICOMB(flagin);<span class="comment">%index corresponding to XIMA and YIMA in the aligned original image vec_A</span>
2770220     vec_B(ind_in)=vec_A(ICOMB);
2780221     vec_B(ind_out)=zeros(size(ind_out));
2790222     A_out{icell}=reshape(vec_B,npy,npx);<span class="comment">%new image in real coordinates</span>
2800223 <span class="keyword">end</span>
2810224 <span class="comment">%Rangx=Rangx-radius_offset;</span>
2820225
2830226</pre></div>
284<hr><address>Generated on Fri 13-Nov-2009 11:17:03 by <strong><a href="http://www.artefact.tk/software/matlab/m2html/">m2html</a></strong> &copy; 2003</address>
285</body>
286</html>
Note: See TracBrowser for help on using the repository browser.