[37] | 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 RUN_STLIN</title> |
---|
| 6 | <meta name="keywords" content="RUN_STLIN"> |
---|
| 7 | <meta name="description" content="'RUN_STLIN': combine velocity fields for stereo PIV"> |
---|
| 8 | <meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1"> |
---|
| 9 | <meta name="generator" content="m2html © 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> > <a href="index.html">.</a> > RUN_STLIN.m</div> |
---|
| 16 | |
---|
| 17 | <!--<table width="100%"><tr><td align="left"><a href="../index.html"><img alt="<" border="0" src="../left.png"> Master index</a></td> |
---|
| 18 | <td align="right"><a href="index.html">Index for . <img alt=">" border="0" src="../right.png"></a></td></tr></table>--> |
---|
| 19 | |
---|
| 20 | <h1>RUN_STLIN |
---|
| 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>'RUN_STLIN': combine velocity fields for stereo PIV</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 RUN_STLIN(file_A,file_B,vel_type,file_st,nx_patch,ny_patch,thresh_patch,fileAxml,fileBxml) </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">'RUN_STLIN': combine velocity fields for stereo PIV |
---|
| 31 | file_A,file_B: input velocity files |
---|
| 32 | vel_type: string ='civ1' or 'civ2'</pre></div> |
---|
| 33 | |
---|
| 34 | <!-- crossreference --> |
---|
| 35 | <h2><a name="_cross"></a>CROSS-REFERENCE INFORMATION <a href="#_top"><img alt="^" border="0" src="../up.png"></a></h2> |
---|
| 36 | This function calls: |
---|
| 37 | <ul style="list-style-image:url(../matlabicon.gif)"> |
---|
| 38 | <li><a href="griddata_uvmat.html" class="code" title="function ZI = griddata_uvmat(X,Y,Z,XI,YI)">griddata_uvmat</a> 'griddata_uvmat': function griddata_uvmat(vec2_X,vec2_Y,vec2_U,vec_X,vec_Y,'linear')</li><li><a href="imadoc2struct.html" class="code" title="function [s,errormsg]=imadoc2struct(ImaDoc)">imadoc2struct</a> 'imadoc2struct': reads the xml file for image documentation</li><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><li><a href="pxcm_tsai.html" class="code" title="function [A11,A12,A13,A21,A22,A23]=pxcm_tsai(a,var_phys)">pxcm_tsai</a> 'pxcm_tsai': find differentials of the Tsai calibration</li><li><a href="read_civxdata.html" class="code" title="function [Field,VelTypeOut]=read_civxdata(filename,FieldNames,VelType)">read_civxdata</a> 'read_civxdata': reads civx data from netcdf files</li><li><a href="struct2nc.html" class="code" title="function errormsg=struct2nc(flname,Data)">struct2nc</a> 'struct2nc': create a netcdf file from a Matlab structure</li><li><a href="warndlg_uvmat.html" class="code" title="function hwarn=warndlg_uvmat(warntext,title)">warndlg_uvmat</a> 'warndlg_uvmat': display warning message (error, warning, confirmation) in a given figure</li></ul> |
---|
| 39 | This function is called by: |
---|
| 40 | <ul style="list-style-image:url(../matlabicon.gif)"> |
---|
| 41 | <li><a href="civ.html" class="code" title="function varargout = civ(varargin)">civ</a> 'civ': function associated with the interface 'civ.fig' for PIV, spline interpolation and stereo PIV (patch)</li></ul> |
---|
| 42 | <!-- crossreference --> |
---|
| 43 | |
---|
| 44 | |
---|
| 45 | <h2><a name="_source"></a>SOURCE CODE <a href="#_top"><img alt="^" border="0" src="../up.png"></a></h2> |
---|
| 46 | <div class="fragment"><pre>0001 <span class="comment">%'RUN_STLIN': combine velocity fields for stereo PIV</span> |
---|
| 47 | 0002 <span class="comment">% file_A,file_B: input velocity files</span> |
---|
| 48 | 0003 <span class="comment">%vel_type: string ='civ1' or 'civ2'</span> |
---|
| 49 | 0004 <a name="_sub0" href="#_subfunctions" class="code">function RUN_STLIN(file_A,file_B,vel_type,file_st,nx_patch,ny_patch,thresh_patch,fileAxml,fileBxml)</a> |
---|
| 50 | 0005 |
---|
| 51 | 0006 [XmlDataA,error]=<a href="imadoc2struct.html" class="code" title="function [s,errormsg]=imadoc2struct(ImaDoc)">imadoc2struct</a>(fileAxml); |
---|
| 52 | 0007 [XmlDataB,error]=<a href="imadoc2struct.html" class="code" title="function [s,errormsg]=imadoc2struct(ImaDoc)">imadoc2struct</a>(fileBxml); |
---|
| 53 | 0008 npxA=[]; npyA=[]; pxB=[]; npyB=[]; |
---|
| 54 | 0009 <span class="keyword">if</span> isfield(XmlDataA,<span class="string">'Camera'</span>) && isfield(XmlDataB,<span class="string">'Camera'</span>) |
---|
| 55 | 0010 <span class="keyword">if</span> isfield(XmlDataA.Camera,<span class="string">'ImageSize'</span>)&& isfield(XmlDataB.Camera,<span class="string">'ImageSize'</span>) |
---|
| 56 | 0011 ImageSizeA=XmlDataA.Camera.ImageSize; |
---|
| 57 | 0012 ImageSizeB=XmlDataB.Camera.ImageSize; |
---|
| 58 | 0013 <span class="keyword">if</span> ~isempty(ImageSizeA)&& ~isempty(ImageSizeB) |
---|
| 59 | 0014 xindex=findstr(ImageSizeA,<span class="string">'x'</span>); |
---|
| 60 | 0015 <span class="keyword">if</span> length(xindex)>=2 |
---|
| 61 | 0016 npxA=str2num(ImageSizeA(1:xindex(1)-1)); |
---|
| 62 | 0017 npyA=str2num(ImageSizeA(xindex(1)+1:xindex(2)-1)); |
---|
| 63 | 0018 <span class="keyword">end</span> |
---|
| 64 | 0019 xindex=findstr(ImageSizeB,<span class="string">'x'</span>); |
---|
| 65 | 0020 <span class="keyword">if</span> length(xindex)>=2 |
---|
| 66 | 0021 npxB=str2num(ImageSizeB(1:xindex(1)-1)); |
---|
| 67 | 0022 npyB=str2num(ImageSizeB(xindex(1)+1:xindex(2)-1)); |
---|
| 68 | 0023 <span class="keyword">end</span> |
---|
| 69 | 0024 <span class="keyword">end</span> |
---|
| 70 | 0025 <span class="keyword">end</span> |
---|
| 71 | 0026 <span class="keyword">end</span> |
---|
| 72 | 0027 <span class="keyword">if</span> isempty(npxA) ||isempty(npxB) |
---|
| 73 | 0028 <a href="warndlg_uvmat.html" class="code" title="function hwarn=warndlg_uvmat(warntext,title)">warndlg_uvmat</a>(<span class="string">'The size of image A needs to be defined in the xml file ImaDoc'</span>,<span class="string">'ERROR'</span>) |
---|
| 74 | 0029 <span class="keyword">return</span> |
---|
| 75 | 0030 <span class="keyword">elseif</span> isempty(npxB) || isempty(npyB) |
---|
| 76 | 0031 <a href="warndlg_uvmat.html" class="code" title="function hwarn=warndlg_uvmat(warntext,title)">warndlg_uvmat</a>(<span class="string">'The size of image B needs to be defined in the xml file ImaDoc'</span>,<span class="string">'ERROR'</span>) |
---|
| 77 | 0032 <span class="keyword">return</span> |
---|
| 78 | 0033 <span class="keyword">end</span> |
---|
| 79 | 0034 <span class="keyword">if</span> isfield(XmlDataA,<span class="string">'GeometryCalib'</span>) |
---|
| 80 | 0035 tsaiA=XmlDataA.GeometryCalib; |
---|
| 81 | 0036 <span class="keyword">else</span> |
---|
| 82 | 0037 <a href="warndlg_uvmat.html" class="code" title="function hwarn=warndlg_uvmat(warntext,title)">warndlg_uvmat</a>(<span class="string">'no geometric calibration available for image A'</span>,<span class="string">'ERROR'</span>) |
---|
| 83 | 0038 <span class="keyword">return</span> |
---|
| 84 | 0039 <span class="keyword">end</span> |
---|
| 85 | 0040 <span class="keyword">if</span> isfield(XmlDataB,<span class="string">'GeometryCalib'</span>) |
---|
| 86 | 0041 tsaiB=XmlDataB.GeometryCalib; |
---|
| 87 | 0042 <span class="keyword">else</span> |
---|
| 88 | 0043 <a href="warndlg_uvmat.html" class="code" title="function hwarn=warndlg_uvmat(warntext,title)">warndlg_uvmat</a>(<span class="string">'no geometric calibration available for image B'</span>,<span class="string">'ERROR'</span>) |
---|
| 89 | 0044 <span class="keyword">return</span> |
---|
| 90 | 0045 <span class="keyword">end</span> |
---|
| 91 | 0046 |
---|
| 92 | 0047 <span class="comment">%corners of each image in real coordinates:</span> |
---|
| 93 | 0048 cornerA(:,1)=[0 0 npxA npxA]';<span class="comment">%x positions</span> |
---|
| 94 | 0049 cornerA(:,2)=[0 npyA 0 npyA]';<span class="comment">%y positions</span> |
---|
| 95 | 0050 cornerB(:,1)=[0 0 npxB npxB]';<span class="comment">%x positions</span> |
---|
| 96 | 0051 cornerB(:,2)=[0 npyB 0 npyB]';<span class="comment">%y positions</span> |
---|
| 97 | 0052 [xyA(:,1),xyA(:,2)]=<a href="phys_XYZ.html" class="code" title="function [Xphys,Yphys,Zphys]=phys_XYZ(Calib,X,Y,Z)">phys_XYZ</a>(tsaiA,cornerA(:,1),cornerA(:,2)); |
---|
| 98 | 0053 [xyB(:,1),xyB(:,2)]=<a href="phys_XYZ.html" class="code" title="function [Xphys,Yphys,Zphys]=phys_XYZ(Calib,X,Y,Z)">phys_XYZ</a>(tsaiB,cornerB(:,1),cornerB(:,2)); |
---|
| 99 | 0054 max_x=max(max(xyA(:,1)),max(xyB(:,1)));<span class="comment">%maximum on the 4 corners of the the images</span> |
---|
| 100 | 0055 min_x=min(min(xyA(:,1)),min(xyB(:,1)));<span class="comment">%minimum on the 4 corners of the the images</span> |
---|
| 101 | 0056 max_y=max(max(xyA(:,2)),max(xyB(:,2))); |
---|
| 102 | 0057 min_y=min(min(xyA(:,2)),min(xyB(:,2))); |
---|
| 103 | 0058 array_realx=[min_x:(max_x-min_x)/(nx_patch-1):max_x]; |
---|
| 104 | 0059 array_realy=[min_y:(max_y-min_y)/(ny_patch-1):max_y]; |
---|
| 105 | 0060 [grid_realx,grid_realy]=meshgrid(array_realx,array_realy); |
---|
| 106 | 0061 grid_real(:,1)=reshape(grid_realx,nx_patch*ny_patch,1); |
---|
| 107 | 0062 grid_real(:,2)=reshape(grid_realy,nx_patch*ny_patch,1); |
---|
| 108 | 0063 grid_real(:,3)=zeros(nx_patch*ny_patch,1); |
---|
| 109 | 0064 [grid_imaA(:,1),grid_imaA(:,2)]=<a href="px_XYZ.html" class="code" title="function [X,Y]=px_XYZ(Calib,Xphys,Yphys,Zphys)">px_XYZ</a>(tsaiA,grid_real(:,1),grid_real(:,2)); |
---|
| 110 | 0065 [grid_imaB(:,1),grid_imaB(:,2)]=<a href="px_XYZ.html" class="code" title="function [X,Y]=px_XYZ(Calib,Xphys,Yphys,Zphys)">px_XYZ</a>(tsaiB,grid_real(:,1),grid_real(:,2)); |
---|
| 111 | 0066 |
---|
| 112 | 0067 flagA=grid_imaA(:,1)>0 & grid_imaA(:,1)<npxA & grid_imaA(:,2)>0 & grid_imaA(:,2)<npyA; |
---|
| 113 | 0068 flagB=grid_imaB(:,1)>0 & grid_imaB(:,1)<npxB & grid_imaB(:,2)>0 & grid_imaB(:,2)<npyB; |
---|
| 114 | 0069 ind_good=find(flagA==1&flagB==1); |
---|
| 115 | 0070 XimaA=grid_imaA(ind_good,1); |
---|
| 116 | 0071 YimaA=grid_imaA(ind_good,2); |
---|
| 117 | 0072 XimaB=grid_imaB(ind_good,1); |
---|
| 118 | 0073 YimaB=grid_imaB(ind_good,2); |
---|
| 119 | 0074 grid_real_x=grid_real(ind_good,1); |
---|
| 120 | 0075 grid_real_y=grid_real(ind_good,2); |
---|
| 121 | 0076 |
---|
| 122 | 0077 <span class="comment">% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span> |
---|
| 123 | 0078 <span class="comment">% %read the velocity fields</span> |
---|
| 124 | 0079 <span class="comment">% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span> |
---|
| 125 | 0080 <span class="comment">%</span> |
---|
| 126 | 0081 <span class="comment">% [dt,time1,pixcmx,pixcmy,vec_X,vec_Y,vec_Z,vec_U,vec_V,vec_W,vec_C,vec_F,fixflag,vel_type_out,error,nb_coord,nb_dim]...</span> |
---|
| 127 | 0082 <span class="comment">% =read_vel({filecell_ncA},{vel_type});</span> |
---|
| 128 | 0083 <span class="comment">%read field A</span> |
---|
| 129 | 0084 [Field,VelTypeOut]=<a href="read_civxdata.html" class="code" title="function [Field,VelTypeOut]=read_civxdata(filename,FieldNames,VelType)">read_civxdata</a>(file_A,[],vel_type); |
---|
| 130 | 0085 <span class="comment">%interpolate on XimaA</span> |
---|
| 131 | 0086 Field.X=Field.X(find(Field.FF==0)); |
---|
| 132 | 0087 Field.Y=Field.Y(find(Field.FF==0)); |
---|
| 133 | 0088 Field.U=Field.U(find(Field.FF==0)); |
---|
| 134 | 0089 Field.V=Field.V(find(Field.FF==0)); |
---|
| 135 | 0090 dXa= <a href="griddata_uvmat.html" class="code" title="function ZI = griddata_uvmat(X,Y,Z,XI,YI)">griddata_uvmat</a>(Field.X,Field.Y,Field.U,XimaA,YimaA); |
---|
| 136 | 0091 dYa= <a href="griddata_uvmat.html" class="code" title="function ZI = griddata_uvmat(X,Y,Z,XI,YI)">griddata_uvmat</a>(Field.X,Field.Y,Field.V,XimaA,YimaA); |
---|
| 137 | 0092 dt=Field.dt; |
---|
| 138 | 0093 time=Field.Time; |
---|
| 139 | 0094 |
---|
| 140 | 0095 <span class="comment">%read field B</span> |
---|
| 141 | 0096 <span class="comment">% [dt,time2,pixcmx,pixcmy,vec_X,vec_Y,vec_Z,vec_U,vec_V,vec_W,vec_C,vec_F,fixflag,vel_type_out,error,nb_coord,nb_dim]...</span> |
---|
| 142 | 0097 <span class="comment">% =read_vel({file_B},{vel_type});</span> |
---|
| 143 | 0098 [Field,VelTypeOut]=<a href="read_civxdata.html" class="code" title="function [Field,VelTypeOut]=read_civxdata(filename,FieldNames,VelType)">read_civxdata</a>(file_B,FieldNames,vel_type); |
---|
| 144 | 0099 <span class="keyword">if</span> ~isequal(Field.dt,dt) |
---|
| 145 | 0100 <a href="warndlg_uvmat.html" class="code" title="function hwarn=warndlg_uvmat(warntext,title)">warndlg_uvmat</a>(<span class="string">'different time intervals for the two velocity fields '</span>,<span class="string">'ERROR'</span>) |
---|
| 146 | 0101 <span class="keyword">return</span> |
---|
| 147 | 0102 <span class="keyword">end</span> |
---|
| 148 | 0103 <span class="keyword">if</span> ~isequal(Field.Time,time) |
---|
| 149 | 0104 <a href="warndlg_uvmat.html" class="code" title="function hwarn=warndlg_uvmat(warntext,title)">warndlg_uvmat</a>(<span class="string">'different times for the two velocity fields '</span>,<span class="string">'ERROR'</span>) |
---|
| 150 | 0105 <span class="keyword">return</span> |
---|
| 151 | 0106 <span class="keyword">end</span> |
---|
| 152 | 0107 <span class="comment">%interpolate on XimaB</span> |
---|
| 153 | 0108 Field.X=Field.X(find(Field.FF==0)); |
---|
| 154 | 0109 Field.Y=Field.Y(find(Field.FF==0)); |
---|
| 155 | 0110 Field.U=Field.U(find(Field.FF==0)); |
---|
| 156 | 0111 Field.V=Field.V(find(Field.FF==0)); |
---|
| 157 | 0112 dXb=<a href="griddata_uvmat.html" class="code" title="function ZI = griddata_uvmat(X,Y,Z,XI,YI)">griddata_uvmat</a>(Field.X,Field.Y,Field.U,XimaB,YimaB); |
---|
| 158 | 0113 dYb=<a href="griddata_uvmat.html" class="code" title="function ZI = griddata_uvmat(X,Y,Z,XI,YI)">griddata_uvmat</a>(Field.X,Field.Y,Field.V,XimaB,YimaB); |
---|
| 159 | 0114 <span class="comment">%eliminate Not-a-Number</span> |
---|
| 160 | 0115 ind_Nan=find(and(~isnan(dXa),~isnan(dXb))); |
---|
| 161 | 0116 dXa=dXa(ind_Nan); |
---|
| 162 | 0117 dYa=dYa(ind_Nan); |
---|
| 163 | 0118 dXb=dXb(ind_Nan); |
---|
| 164 | 0119 dYb=dYb(ind_Nan); |
---|
| 165 | 0120 grid_phys1(:,1)=grid_real_x(ind_Nan); |
---|
| 166 | 0121 grid_phys1(:,2)=grid_real_y(ind_Nan); |
---|
| 167 | 0122 |
---|
| 168 | 0123 <span class="comment">%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span> |
---|
| 169 | 0124 <span class="comment">%compute the coefficients</span> |
---|
| 170 | 0125 <span class="comment">%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span> |
---|
| 171 | 0126 |
---|
| 172 | 0127 |
---|
| 173 | 0128 [A11,A12,A13,A21,A22,A23]=<a href="pxcm_tsai.html" class="code" title="function [A11,A12,A13,A21,A22,A23]=pxcm_tsai(a,var_phys)">pxcm_tsai</a>(tsaiA,grid_phys1); |
---|
| 174 | 0129 [B11,B12,B13,B21,B22,B23]=<a href="pxcm_tsai.html" class="code" title="function [A11,A12,A13,A21,A22,A23]=pxcm_tsai(a,var_phys)">pxcm_tsai</a>(tsaiB,grid_phys1); |
---|
| 175 | 0130 |
---|
| 176 | 0131 C1=A11.*A22-A12.*A21; |
---|
| 177 | 0132 C2=A13.*A22-A12.*A23; |
---|
| 178 | 0133 C3=A13.*A21-A11.*A23; |
---|
| 179 | 0134 D1=B11.*B22-B12.*B21; |
---|
| 180 | 0135 D2=B13.*B22-B12.*B23; |
---|
| 181 | 0136 D3=B13.*B21-B11.*B23; |
---|
| 182 | 0137 A1=(A22.*D1.*(C1.*D3-C3.*D1)+A21.*D1.*(C2.*D1-C1.*D2)); |
---|
| 183 | 0138 A2=(A12.*D1.*(C3.*D1-C1.*D3)+A11.*D1.*(C1.*D2-C2.*D1)); |
---|
| 184 | 0139 B1=(B22.*C1.*(C3.*D1-C1.*D3)+B21.*C1.*(C1.*D2-C2.*D1)); |
---|
| 185 | 0140 B2=(B12.*C1.*(C1.*D3-C3.*D1)+B11.*C1.*(C2.*D1-C1.*D2)); |
---|
| 186 | 0141 Lambda=(A1.*dXa+A2.*dYa+B1.*dXb+B2.*dYb)./(A1.*A1+A2.*A2+B1.*B1+B2.*B2); |
---|
| 187 | 0142 |
---|
| 188 | 0143 <span class="comment">%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span> |
---|
| 189 | 0144 <span class="comment">%Projection for compatible displacements</span> |
---|
| 190 | 0145 <span class="comment">%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span> |
---|
| 191 | 0146 Ua=dXa-Lambda.*A1; |
---|
| 192 | 0147 Va=dYa-Lambda.*A2; |
---|
| 193 | 0148 Ub=dXb-Lambda.*B1; |
---|
| 194 | 0149 Vb=dYb-Lambda.*B2; |
---|
| 195 | 0150 |
---|
| 196 | 0151 <span class="comment">%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span> |
---|
| 197 | 0152 <span class="comment">%Calculations of displacements and error</span> |
---|
| 198 | 0153 <span class="comment">%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%</span> |
---|
| 199 | 0154 U=(A22.*D2.*Ua-A12.*D2.*Va-B22.*C2.*Ub+B12.*C2.*Vb)./(C1.*D2-C2.*D1); |
---|
| 200 | 0155 V=(A21.*D3.*Ua-A11.*D3.*Va-B21.*C3.*Ub+B11.*C3.*Vb)./(C3.*D1-C1.*D3); |
---|
| 201 | 0156 W=(A22.*D1.*Ua-A12.*D1.*Va-B22.*C1.*Ub+B12.*C1.*Vb)./(C2.*D1-C1.*D2); |
---|
| 202 | 0157 W1=(-A21.*D1.*Ua+A11.*D1.*Va+B21.*C1.*Ub-B11.*C1.*Vb)./(C1.*D3-C3.*D1); |
---|
| 203 | 0158 |
---|
| 204 | 0159 error=sqrt((A1.*dXa+A2.*dYa+B1.*dXb+B2.*dYb).*(A1.*dXa+A2.*dYa+B1.*dXb+B2.*dYb)./(A1.*A1+A2.*A2+B1.*B1+B2.*B2)); |
---|
| 205 | 0160 |
---|
| 206 | 0161 ind_error=(find(error<thresh_patch)); |
---|
| 207 | 0162 U=U(ind_error); |
---|
| 208 | 0163 V=V(ind_error); |
---|
| 209 | 0164 W=W(ind_error);<span class="comment">%correction for water interface</span> |
---|
| 210 | 0165 error=error(ind_error); |
---|
| 211 | 0166 |
---|
| 212 | 0167 <span class="comment">%create nc grid file</span> |
---|
| 213 | 0168 Result.ListGlobalAttribute={<span class="string">'nb_coord'</span>,<span class="string">'nb_dim'</span>,<span class="string">'constant_pixcm'</span>,<span class="string">'absolut_time_T0'</span>,<span class="string">'hart'</span>,<span class="string">'dt'</span>,<span class="string">'civ'</span>}; |
---|
| 214 | 0169 Result.nb_coord=3;<span class="comment">%grid file, no velocity</span> |
---|
| 215 | 0170 Result.nb_dim=2; |
---|
| 216 | 0171 Result.constant_pixcm=0;<span class="comment">%no linear correspondance with images</span> |
---|
| 217 | 0172 Result.absolut_time_T0=time;<span class="comment">%absolute time of the field</span> |
---|
| 218 | 0173 Result.hart=0; |
---|
| 219 | 0174 Result.dt=dt;<span class="comment">%time interval for image correlation (put by default)</span> |
---|
| 220 | 0175 <span class="comment">% cte.title='grid';</span> |
---|
| 221 | 0176 Result.civ=0;<span class="comment">%not a civ file (no direct correspondance with an image)</span> |
---|
| 222 | 0177 Result.ListDimName={<span class="string">'nb_vectors'</span>} |
---|
| 223 | 0178 Result.DimValue=length(U); |
---|
| 224 | 0179 Result.ListVarName={<span class="string">'vec_X'</span>;<span class="string">'vec_Y'</span>;<span class="string">'vec_U'</span>;<span class="string">'vec_V'</span>;<span class="string">'vec_W'</span>;<span class="string">'vec_E'</span>}; |
---|
| 225 | 0180 Result.VarDimIndex: {[1] [1] [1] [1] [1] [1]} |
---|
| 226 | 0181 Result.vec_X= grid_phys1(ind_error,1); |
---|
| 227 | 0182 Result.vec_Y= grid_phys1(ind_error,2); |
---|
| 228 | 0183 Result.vec_U=U/dt; |
---|
| 229 | 0184 Result.vec_V=V/dt; |
---|
| 230 | 0185 Result.vec_W=W/dt; |
---|
| 231 | 0186 Result.vec_E=error; |
---|
| 232 | 0187 <span class="comment">% error=write_netcdf(file_st,cte,fieldlabels,grid_phys);</span> |
---|
| 233 | 0188 error=<a href="struct2nc.html" class="code" title="function errormsg=struct2nc(flname,Data)">struct2nc</a>(file_st,Result); |
---|
| 234 | 0189 display([file_st <span class="string">' written'</span>]) |
---|
| 235 | 0190 |
---|
| 236 | 0191 |
---|
| 237 | 0192 |
---|
| 238 | 0193 |
---|
| 239 | 0194 |
---|
| 240 | 0195 |
---|
| 241 | 0196 |
---|
| 242 | 0197 |
---|
| 243 | 0198 |
---|
| 244 | 0199 |
---|
| 245 | 0200 |
---|
| 246 | 0201 |
---|
| 247 | 0202</pre></div> |
---|
| 248 | <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> © 2003</address> |
---|
| 249 | </body> |
---|
| 250 | </html> |
---|