function beam_explorer( h_arrayant, Jp )
%BEAM_EXPLORER Creates an interactive plot of the beam-forming capabilities of an array antenna
%
% Calling object:
%   Single object
%
% Description:
%   When applying maximum-ratio transmission (MRT) to calculate the coupling weights of an array
%   antenna, it is possible to direct a beam towards a given direction. However, the antenna
%   geometry (i.e. the positions and the orientations of the individual elements) and the shape of
%   the individual element-patterns will determine the overall shape of the beam and the existence
%   and magnitude of so-called sidelobes.  This method creates an interactive plot that uses the
%   mouse pointer position to determine the target direction. The y-axis corresponds to the
%   elevation direction and the x-axis corresponds to the azimuth direction relative to the local
%   antenna coordinate system. Then, the method calculates the MRT-weights that direct the beam
%   towards this position and applies it to the antenna pattern. The plot is then updated in real-
%   time to visualize the radiated power using MRT beamforming. The maximum is normalized to 1, the
%   minimum is normalized to 0. Ideally, the array antenna creates a single narrow beam that
%   coincides exactly with the target direction, i.e. the maximum is always under the mouse
%   pointer. However, design limitations (i.e., using planar or circular arrays, number of
%   elements, etc.) will either lead to unwanted side-lobes or a widening of the main lobe. By
%   clicking the left mouse button, the animation is paused and it is possible to use the data-
%   pointer to read the values from the plot.  
%
%   Note: It is important to set the correct center frequency in the array object.
%
% Input:
%   Jp
%   The polarization (Jones-vector) of the probe antenna, Default [ 1 ; 0 ]
%
%
% QuaDRiGa Copyright (C) 2011-2019
% Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V. acting on behalf of its
% Fraunhofer Heinrich Hertz Institute, Einsteinufer 37, 10587 Berlin, Germany
% All rights reserved.
%
% e-mail: quadriga@hhi.fraunhofer.de
%
% This file is part of QuaDRiGa.
%
% The Quadriga software is provided by Fraunhofer on behalf of the copyright holders and
% contributors "AS IS" and WITHOUT ANY EXPRESS OR IMPLIED WARRANTIES, including but not limited to
% the implied warranties of merchantability and fitness for a particular purpose.
%
% You can redistribute it and/or modify QuaDRiGa under the terms of the Software License for 
% The QuaDRiGa Channel Model. You should have received a copy of the Software License for The
% QuaDRiGa Channel Model along with QuaDRiGa. If not, see <http://quadriga-channel-model.de/>. 

global b Wa Wb im

if ~exist('Jp','var') || isempty( Jp )
   Jp = [1;0];
end

if numel( h_arrayant ) > 1
   error('QuaDRiGa:qd_arrayant:beam_explorer','beam_explorer is not definded for object arrays.');
else
    h_arrayant = h_arrayant(1,1); % workaround for octave
end

tmp = h_arrayant.copy;
tmp.combine_pattern;

Wa = single( reshape( tmp.Fa, [], h_arrayant.no_elements ) );
Wb = single( reshape( tmp.Fb, [], h_arrayant.no_elements ) );

b = qd_builder('LOSonly');
b.simpar.center_frequency = h_arrayant.center_frequency;
b.simpar.show_progress_bars = 0;
b.tx_array = h_arrayant;
b.rx_array = qd_arrayant('omni');
b.rx_array.Fa(:) = Jp(1);
b.rx_array.Fb(:) = Jp(2);
b.tx_position = [0;0;0];
b.rx_positions = [100;0;0];
b.gen_parameters;

if b.NumClusters == 1
    c = b.get_los_channels('single','coeff');
else
    c = b.get_channels;
    c = single( sum(c.coeff,3) );
end
cf = permute( c, [2,3,1] ) ;

A = abs( Wa * conj(cf) ).^2 + abs( Wb * conj(cf) ).^2;
A = A - min(A);
A = A./max(A);
A = reshape( A , h_arrayant.no_el , h_arrayant.no_az , [] );

im = imagesc(h_arrayant.azimuth_grid*180/pi,h_arrayant.elevation_grid*180/pi,A);
set(gca,'Ydir','Normal');
set(gca,'clim',[0 1]);
set(gca,'xlim',[-180 180]);
set(gca,'ylim',[-90 90]);
set(gcf,'Pointer','crosshair')
set(gcf,'WindowButtonMotionFcn', @mouseMove);
set(gcf,'WindowButtonDownFcn', @mouseClick);

end

function mouseMove(object, eventdata)

global b Wa Wb im 

C = get(gca, 'CurrentPoint');
az = C(1,1);
el = C(1,2);

b.rx_positions = [ cosd(az)*cosd(el) ; sind(az)*cosd(el) ; sind(el) ]*100;
c = b.get_los_channels('single','coeff');
cf = permute( c, [2,3,1] ) ;

A = abs( Wa * conj(cf) ).^2 + abs( Wb * conj(cf) ).^2;
A = A - min(A);
A = A./max(A);
A = reshape( A , b.tx_array.no_el , b.tx_array.no_az , [] );

set( im, 'CData', A);

title(gca, ['(az,el) = (', num2str(az), ', ',num2str(el), ')']);

end

function mouseClick(object, eventdata)

persistent stp
if isempty( stp )
   stp = false;
end

if stp
    set(gcf,'WindowButtonMotionFcn', @mouseMove);
    stp = ~stp;
else
    set(gcf,'WindowButtonMotionFcn', '');
    stp = ~stp;
end

end
