source: trunk/src/proj_grid.m @ 965

Last change on this file since 965 was 924, checked in by g7moreau, 8 years ago
  • Update Copyright Copyright 2008-2016, LEGI UMR 5519 / CNRS UGA G-INP, Grenoble, France
File size: 5.1 KB
Line 
1%'proj_grid': project  fields with unstructured coordinantes on a regular grid
2% -------------------------------------------------------------------------
3% function [A,rangx,rangy]=proj_grid(vec_X,vec_Y,vec_A,rgx_in,rgy_in,npxy_in)
4
5%=======================================================================
6% Copyright 2008-2016, LEGI UMR 5519 / CNRS UGA G-INP, Grenoble, France
7%   http://www.legi.grenoble-inp.fr
8%   Joel.Sommeria - Joel.Sommeria (A) legi.cnrs.fr
9%
10%     This file is part of the toolbox UVMAT.
11%
12%     UVMAT is free software; you can redistribute it and/or modify
13%     it under the terms of the GNU General Public License as published
14%     by the Free Software Foundation; either version 2 of the license,
15%     or (at your option) any later version.
16%
17%     UVMAT is distributed in the hope that it will be useful,
18%     but WITHOUT ANY WARRANTY; without even the implied warranty of
19%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
20%     GNU General Public License (see LICENSE.txt) for more details.
21%=======================================================================
22
23function [A,rangx,rangy]=proj_grid(vec_X,vec_Y,vec_A,rgx_in,rgy_in,npxy_in)
24    if length(vec_Y)<2
25        msgbox_uvmat('ERROR','less than 2 points in proj_grid.m');
26        return;
27    end
28    diffy=diff(vec_Y); %difference dy=vec_Y(i+1)-vec_Y(i)
29    index=find(diffy);% find the indices of vec_Y after wich a change of horizontal line occurs(diffy non zero)
30    if isempty(index); msgbox_uvmat('ERROR','points aligned along abscissa in proj_grid.m'); return; end;%points aligned% A FAIRE: switch to line plot.
31    diff2=diff(diffy(index));% diff2 = fluctuations of the detected vertical grid mesh dy
32    if max(abs(diff2))>0.001*abs(diffy(index(1))) % if max(diff2) is larger than 1/1000 of the first mesh dy
33        % the data are not regularly spaced and must be interpolated  on a regular grid
34        if exist('rgx_in','var') & ~isempty (rgx_in) & isnumeric(rgx_in) & length(rgx_in)==2%  positions imposed from input
35            rangx=rgx_in; % first and last positions
36            rangy=rgy_in;
37%             npxy=npxy_in;
38            dxy(1)=1/(npxy_in(1)-1);%grid mesh in y
39            dxy(2)=1/(npxy_in(2)-1);%grid mesh in x
40            dxy(1)=(rangy(2)-rangy(1))/(npxy_in(1)-1);%grid mesh in y
41            dxy(2)=(rangx(2)-rangx(1))/(npxy_in(2)-1);%grid mesh in x
42        else % interpolation grid automatically determined
43            rangx(1)=min(vec_X);
44            rangx(2)=max(vec_X);
45            rangy(2)=min(vec_Y);
46            rangy(1)=max(vec_Y);
47            dxymod=sqrt((rangx(2)-rangx(1))*(rangy(1)-rangy(2))/length(vec_X));
48            dxy=[-dxymod/4 dxymod/4];% increase the resolution 4 times
49        end
50        xi=[rangx(1):dxy(2):rangx(2)];
51        yi=[rangy(1):dxy(1):rangy(2)];
52        [XI,YI]=meshgrid(xi,yi);% creates the matrix of regular coordinates
53        A=griddata_uvmat(vec_X,vec_Y,vec_A,xi,yi');
54        A=reshape(A,length(yi),length(xi));
55    else
56        x=vec_X(1:index(1));% the set of abscissa (obtained on the first line)
57        indexend=index(end);% last vector index of line change
58        ymax=vec_Y(indexend+1);% y coordinate AFTER line change
59        ymin=vec_Y(index(1));
60        %y=[vec_Y(index) ymax]; % the set of y ordinates including the last one
61        y=vec_Y(index);
62        y(length(y)+1)=ymax;
63        nx=length(x);   %number of grid points in x
64        ny=length(y);   % number of grid points in y
65        B=(reshape(vec_A,nx,ny))'; %vec_A reshaped as a rectangular matrix
66        [X,Y]=meshgrid(x,y);% positions X and Y also reshaped as matrix
67       
68        %linear interpolation to improve the image resolution and/or adjust
69        %to prescribed positions
70        test_interp=1;
71        if exist('rgx_in','var') & ~isempty (rgx_in) & isnumeric(rgx_in) & length(rgx_in)==2%  positions imposed from input
72            rangx=rgx_in; % first and last positions
73            rangy=rgy_in;
74            npxy=npxy_in;
75        else       
76            rangx=[vec_X(1) vec_X(nx)];% first and last position found for x
77%             rangy=[ymin ymax];
78              rangy=[max(ymax,ymin) min(ymax,ymin)];
79            if max(nx,ny) <= 64 & isequal(npxy_in,'np>256')
80                npxy=[8*ny 8*nx];% increase the resolution 8 times
81            elseif max(nx,ny) <= 128 & isequal(npxy_in,'np>256')
82                npxy=[4*ny 4*nx];% increase the resolution 4 times
83            elseif max(nx,ny) <= 256 & isequal(npxy_in,'np>256')
84                npxy=[2*ny 2*nx];% increase the resolution 2 times
85            else
86                npxy=[ny nx];
87                test_interp=0; % no interpolation done
88            end
89        end
90        if test_interp==1%if we interpolate
91            xi=[rangx(1):(rangx(2)-rangx(1))/(npxy(2)-1):rangx(2)];
92            yi=[rangy(1):(rangy(2)-rangy(1))/(npxy(1)-1):rangy(2)];
93            [XI,YI]=meshgrid(xi,yi);
94            A = interp2(X,Y,B,XI,YI);
95        else %no interpolation for a resolution higher than 256
96            A=B;
97            XI=X;
98            YI=Y;
99        end
100    end
Note: See TracBrowser for help on using the repository browser.