ISOBSERVABLE Test for observability. returns 1 if the system is observable, 0 otherwise. created by Rui Bo on Jan 9, 2010
0001 function TorF = isobservable(H, pv, pq) 0002 %ISOBSERVABLE Test for observability. 0003 % returns 1 if the system is observable, 0 otherwise. 0004 % created by Rui Bo on Jan 9, 2010 0005 0006 % MATPOWER 0007 % $Id: isobservable.m,v 1.3 2010/04/26 19:45:26 ray Exp $ 0008 % by Rui Bo 0009 % Copyright (c) 2009-2010 by Rui Bo 0010 % 0011 % This file is part of MATPOWER. 0012 % See http://www.pserc.cornell.edu/matpower/ for more info. 0013 % 0014 % MATPOWER is free software: you can redistribute it and/or modify 0015 % it under the terms of the GNU General Public License as published 0016 % by the Free Software Foundation, either version 3 of the License, 0017 % or (at your option) any later version. 0018 % 0019 % MATPOWER is distributed in the hope that it will be useful, 0020 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0021 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0022 % GNU General Public License for more details. 0023 % 0024 % You should have received a copy of the GNU General Public License 0025 % along with MATPOWER. If not, see <http://www.gnu.org/licenses/>. 0026 % 0027 % Additional permission under GNU GPL version 3 section 7 0028 % 0029 % If you modify MATPOWER, or any covered work, to interface with 0030 % other modules (such as MATLAB code and MEX-files) available in a 0031 % MATLAB(R) or comparable environment containing parts covered 0032 % under other licensing terms, the licensors of MATPOWER grant 0033 % you additional permission to convey the resulting work. 0034 0035 %% options 0036 tol = 1e-5; % mpopt(2); 0037 check_reason = 1; % check reason for system being not observable 0038 % 0: no check 0039 % 1: check (NOTE: may be time consuming due to svd calculation) 0040 0041 %% test if H is full rank 0042 [m, n] = size(H); 0043 r = rank(H); 0044 if r < min(m, n) 0045 TorF = 0; 0046 else 0047 TorF = 1; 0048 end 0049 0050 %% look for reasons for system being not observable 0051 if check_reason && ~TorF 0052 %% look for variables not being observed 0053 idx_trivialColumns = []; 0054 varNames = {}; 0055 for j = 1:n 0056 normJ = norm(H(:, j), inf); 0057 if normJ < tol % found a zero column 0058 idx_trivialColumns = [idx_trivialColumns j]; 0059 varName = getVarName(j, pv, pq); 0060 varNames{length(idx_trivialColumns)} = varName; 0061 end 0062 end 0063 0064 if ~isempty(idx_trivialColumns) % found zero-valued column vector 0065 fprintf('Warning: The following variables are not observable since they are not related with any measurement!'); 0066 varNames 0067 idx_trivialColumns 0068 else % no zero-valued column vector 0069 %% look for dependent column vectors 0070 for j = 1:n 0071 rr = rank(H(:, 1:j)); 0072 if rr ~= j % found dependent column vector 0073 %% look for linearly depedent vector 0074 colJ = H(:, j); % j(th) column of H 0075 varJName = getVarName(j, pv, pq); 0076 for k = 1:j-1 0077 colK = H(:, k); 0078 if rank([colK colJ]) < 2 % k(th) column vector is linearly dependent of j(th) column vector 0079 varKName = getVarName(k, pv, pq); 0080 fprintf('Warning: %d(th) column vector (w.r.t. %s) of H is linearly dependent of %d(th) column vector (w.r.t. %s)!\n', j, varJName, k, varKName); 0081 return; 0082 end 0083 end 0084 end 0085 end 0086 fprintf('Warning: No specific reason was found for system being not observable.\n'); 0087 end 0088 end