-
Notifications
You must be signed in to change notification settings - Fork 0
/
reale.m
98 lines (87 loc) · 1.95 KB
/
reale.m
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
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
function y = reale(x, arg2, arg3)
%|
%| y = reale(x)
%| y = reale(x, tol)
%| y = reale(x, 'warn', 'message')
%| y = reale(x, 'error')
%| y = reale(x, 'report')
%| y = reale(x, 'prompt')
%| y = reale(x, 'disp')
%|
%| return real part of complex data (with error checking).
%| checks that imaginary part is negligible (or warning etc. if not)
%|
%| Copyright Jeff Fessler, University of Michigan
if nargin < 1, help(mfilename), error(mfilename), end
if nargin == 1 && streq(x, 'test'), reale_test, return, end
com = 'error';
if isa(x, 'double')
tol = 1e-13;
else
tol = 1e-6;
end
if nargin > 1
if ischar(arg2)
com = arg2;
elseif isnumeric(arg2)
tol = arg2;
end
end
if streq(com, 'disp')
;
elseif streq(com, 'warn')
onlywarn = 1;
elseif streq(com, 'error')
onlywarn = 0;
elseif streq(com, 'prompt')
;
elseif streq(com, 'report')
;
else
fail('bad argument %s', com)
end
max_abs_x = max(abs(x(:)));
if max_abs_x == 0
if any(imag(x(:)) ~= 0)
error 'max real 0, but imaginary!'
else
y = real(x);
return
end
end
frac = max(abs(imag(x(:)))) / max_abs_x;
if streq(com, 'report')
printm('imaginary part %g%%', frac * 100)
return
end
if frac > tol
[cname line] = caller_name;
t = sprintf('%s(%d): %s: imaginary fraction of %s [class %s] is %g', ...
cname, line, mfilename, inputname(1), class(x), frac);
if isvar('arg3')
t = [t ', ' arg3];
end
if streq(com, 'disp')
disp(t)
elseif streq(com, 'prompt')
printm('reale() called for input with imaginary part %g%%', frac * 100)
printm('reale() called in context where a large imaginary part')
printm('is likely an *error*. proceed with caution!')
t = input('proceed? [y|n]: ', 's');
if isempty(t) || t(1) ~= 'y'
printm('ok, aborting is probably wise!')
error ' '
end
elseif onlywarn
disp(t)
else
error(t)
end
end
y = real(x);
function reale_test
x = 7 + 1i*eps;
reale(x, 'warn');
reale(x, 'prompt');
%reale(x, 'report'); % check error reporting
%reale(x, eps/100) % check error faulting