function [map,ens] = eigmode_map(wf,eig_en,en_res)
% eigmode_map  - construct eigenmode visualization map
%   Usage:
% Input:
%   wf - array of eigenmode wave functions of the size length(eig_en)x(system length)
%   eig_en - eigenmode energies
%   en_res - energy resolution for the visualization
% Output:
%   map - visualization matrix of the size Nx(system length), where
%       N=5*round((max(eig_en)-min(eig_en))/en_res)
%   ens - vector of energies of the rows in 
   
% map energies and # of rows in m:
subres=2; %5;
del=0.2; % Delta energy in units of eig_en
ens=[(min(eig_en)-del):(en_res/subres):max(eig_en)]; % Energies of the map
N=length(ens);

% Weight matrix (column is map energy, row is eigen-energy)
w=exp(-(repmat(ens',1,length(eig_en))-repmat(eig_en,N,1)).^2/en_res^2);

map=w*(abs(wf)).^2; % Some intensities (incoherently)
% map=abs(w*wf).^2; % Sum amplitudes (coherently)