dworkLib/isHermitian.m at master · LarsonLab/dworkLib · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61

function [out,error] = isHermitian( data, varargin )
  % [out,error] = isHermitian( data, [ threshold ] )
  % Determines whether or not the real part of data is even and the
  %   imaginary part is odd (with circular boundary conditions)
  % Indexes are defined according to fftshift
  %
  % Inputs:
  % data is a 1D or 2D array
  % threshold is an optional input.  The relative error
  %   must be less than threshold to be considered Hermitian.

  defaultThresh = 0;
  p = inputParser;
  p.addOptional( 'thresh', defaultThresh, @isnumeric );
  p.parse( varargin{:} );
  thresh = p.Results.thresh;

  numDims = ndims( data );

  if isrow(data) || iscolumn(data)
    [out,error] = isHermitian1D( data, thresh );

  elseif numDims == 2
    [out,error] = isHermitian2D( data, thresh );

  end

end


function [out,error] = isHermitian1D( data, thresh )
  mirrorData = flipud( data(:) );

  nData = numel( data );

  nEven = ~mod(nData,2);
  if nEven
    mirrorData = circshift( mirrorData, [nEven 1] );
  end

  diff = data - conj(mirrorData);
  error = norm(diff(:),2) / norm(data(:),2);
  if error > thresh, out = false; else out = true; end;
end


function [out,error] = isHermitian2D( img, thresh )
  mirrorImg = rot90( img, 2 );

  [ny,nx] = size(img);
  nyEven = ~mod(ny,2);
  nxEven = ~mod(nx,2);
  if nyEven || nxEven
    mirrorImg = circshift( mirrorImg, +[nyEven nxEven] );
  end

  diff = img - conj(mirrorImg);
  error = norm(diff(:),2) / norm(img(:),2);
  if error > thresh, out = false; else out = true; end;
end