| [510] | 1 | % 'FFT': calculate and display 2D spectrum of the input scalar
|
|---|
| 2 | % GUI_input=FFT(hget_field)
|
|---|
| 3 | %
|
|---|
| 4 | % OUTPUT:
|
|---|
| 5 | % GUI_input: option for display in the GUI get_field
|
|---|
| 6 | %
|
|---|
| 7 | %INPUT:
|
|---|
| 8 | % hget_field: handles of the GUI get_field
|
|---|
| 9 | %
|
|---|
| 10 |
|
|---|
| 11 | function DataOut=FFT2_detrend(DataIn)
|
|---|
| 12 | %% set GUI config
|
|---|
| 13 | DataOut=[];
|
|---|
| 14 | if strcmp(DataIn,'*')
|
|---|
| 15 | DataOut.InputFieldType='scalar';
|
|---|
| 16 | return
|
|---|
| 17 | end
|
|---|
| 18 | %%%%%%%%%%%%%%%%%%%%%%%%%
|
|---|
| [513] | 19 | [CellVarIndex,NbDim,CellVarType,errormsg]=find_field_cells(DataIn);
|
|---|
| [510] | 20 | if ~isempty(errormsg)
|
|---|
| 21 | DataOut.Txt=errormsg;
|
|---|
| 22 | return
|
|---|
| 23 | end
|
|---|
| 24 | DataOut.ListVarName={};
|
|---|
| 25 | DataOut.VarDimName={};
|
|---|
| 26 | for ilist=1:numel(CellVarIndex)
|
|---|
| 27 | if NbDim(ilist)==2 && numel(CellVarType{ilist}.coord)==2 % field with structured coordinates
|
|---|
| 28 | %process coordinates
|
|---|
| 29 | CoordName=DataIn.ListVarName(CellVarType{ilist}.coord);
|
|---|
| 30 | x1 = DataIn.(CoordName{2}); y1 = DataIn.(CoordName{1});
|
|---|
| 31 | [x y] = meshgrid(x1,y1);
|
|---|
| 32 | coeff(1,1) = sum(sum(x.^2)); coeff(1,2) = sum(sum(x.*y)); coeff(1,3) = sum(sum(x));
|
|---|
| 33 | coeff(2,1) = sum(sum(x.*y)); coeff(2,2) = sum(sum(y.^2)); coeff(2,3) = sum(sum(y));
|
|---|
| 34 | coeff(3,1) = sum(sum(x)); coeff(3,2) = sum(sum(y)); coeff(3,3) = length(x1)*length(y1);
|
|---|
| 35 | delta_x = x1(2) - x1(1); delta_y = y1(2) - y1(1);
|
|---|
| 36 | Nx = length(x1); Ny = length(y1);
|
|---|
| 37 | Nxa = 1:Nx; Nya = 1:Ny;
|
|---|
| 38 | ssx = find(Nxa<Nx/2); ssy = find(Nya<Ny/2);
|
|---|
| 39 | Nxa(Nx-ssx+1) = -Nxa(ssx)+1; Nya(Ny-ssy+1) = -Nya(ssy)+1;
|
|---|
| 40 | [Nxa Ix] = sort(Nxa); [Nya Iy] = sort(Nya);
|
|---|
| 41 | kx1 = (2*pi/delta_x/Nx)*(Nxa-1); ky1 = (2*pi/delta_y/Ny)*(Nya-1);
|
|---|
| 42 | ss = find(ky1>=0); ky1 = ky1(ss);
|
|---|
| 43 | [kx ky] = meshgrid(kx1,ky1);
|
|---|
| 44 | DataOut.(CoordName{2}) = kx1; DataOut.(CoordName{1}) = ky1;
|
|---|
| 45 | if isfield(DataIn,'CoordUnit')
|
|---|
| 46 | DataOut.CoordUnit=[DataIn.CoordUnit '^{-1}'];
|
|---|
| 47 | DataOut.ListGlobalAttribute={'CoordUnit'};
|
|---|
| 48 | end
|
|---|
| 49 | %process scalar
|
|---|
| 50 | ivar=CellVarType{ilist}.scalar(1);
|
|---|
| 51 | VarName=DataIn.ListVarName{ivar};
|
|---|
| 52 | z=DataIn.(VarName);
|
|---|
| 53 | rhs(1) = sum(sum(x.*z)); rhs(2) = sum(sum(y.*z)); rhs(3) = sum(sum(z));
|
|---|
| 54 | lin_coeff = inv(coeff)*rhs';
|
|---|
| 55 | lin_trend = lin_coeff(1)*x + lin_coeff(2)*y + lin_coeff(3);
|
|---|
| 56 | z2 = z - lin_trend;
|
|---|
| 57 | spec2 = abs(fft2(z2)).^2;
|
|---|
| 58 | spec2 = spec2(Iy,Ix);
|
|---|
| 59 | spec2 = spec2(ss,:);
|
|---|
| 60 | %DataOut.(VarName) = log(spec2);
|
|---|
| 61 | DataOut.(VarName) = spec2;
|
|---|
| 62 | spec_sum=sum(sum(spec2));
|
|---|
| 63 | kx_mean=sum(sum(spec2.*kx))/spec_sum
|
|---|
| 64 | ky_mean=sum(sum(spec2.*ky))/spec_sum
|
|---|
| 65 | theta=atand(ky_mean/kx_mean)
|
|---|
| 66 | lambda=2*pi/(sqrt(kx_mean*kx_mean+ky_mean*ky_mean))
|
|---|
| 67 | %DataOut.ListVarName=[CoordName {VarName} {'kx'} {'ky'}];%list of variables
|
|---|
| 68 | %DataOut.VarDimName=[CoordName {CoordName} {'one'} {'one'}];%list of dimensions for variables
|
|---|
| 69 |
|
|---|
| 70 | DataOut.ListVarName=[CoordName {VarName}];%list of variables
|
|---|
| 71 | DataOut.VarDimName=[CoordName {CoordName}];%list of dimensions for variables
|
|---|
| 72 | break
|
|---|
| 73 | end
|
|---|
| 74 | end
|
|---|
| 75 |
|
|---|
| 76 |
|
|---|
| 77 |
|
|---|
