Search code examples
matlabplotsignal-processingmatlab-figuremodulation

How to plot colorful histogram type constellation diagram in Matlab


I would like to plot constellation diagram similar to the figure below. Colorful Constellation.

My approach is something like this

 clc;
 clear all;
 close all;
 N=30000;                            
 M=16;                               
 Sr=randint(N,1,[0,(M-1)]);          
 S=qammod(Sr,16,0,'gray'); S=S(:);   
 Noisy_Data=awgn(S,20,'measured');       % Add AWGN
 figure(2)
 subplot(1,2,1)
 plot(S,'o','markersize',10);
 grid on
 subplot(1,2,2)
 plot(Noisy_Data,'.');
 grid on

May you assist me to make necessary modification to get graph similar to the figure attached above. Thank you.


Solution

  • The first thing to do would be to compute a 2D histogram of your data. This can be done with the following:

    % Size of the histogram matrix
    Nx   = 160;
    Ny   = 160;
    
    % Choose the bounds of the histogram to match min/max of data samples.
    % (you could alternatively use fixed bound, e.g. +/- 4)
    ValMaxX = max(real(Noisy_Data));
    ValMinX = min(real(Noisy_Data));
    ValMaxY = max(imag(Noisy_Data));
    ValMinY = min(imag(Noisy_Data));
    dX = (ValMaxX-ValMinX)/(Nx-1);
    dY = (ValMaxY-ValMinY)/(Ny-1);
    
    % Figure out which bin each data sample fall into
    IdxX = 1+floor((real(Noisy_Data)-ValMinX)/dX);
    IdxY = 1+floor((imag(Noisy_Data)-ValMinY)/dY);
    H = zeros(Ny,Nx);
    for i=1:N
      if (IdxX(i) >= 1 && IdxX(i) <= Nx && IdxY(i) >= 1 && IdxY(i) <= Ny)
        % Increment histogram count
        H(IdxY(i),IdxX(i)) = H(IdxY(i),IdxX(i)) + 1;
      end
    end
    

    Note that you can play around with parameters Nx and Ny to adjust the desired resolution of the plot. Keep in mind that the larger the histogram, the more data samples (controlled by the parameter N of your simulation) you'll need to have enough data in the histogram bins to avoid getting a spotty plot.

    You can then plot the histogram as a color map based on this answer. In doing so, you likely would want to add a constant to all non-zero bins of the histogram so that the white band is reserved for zero valued bins. This would provide a better correlation with the scatter plot. This can be done with:

    % Colormap that approximate the sample figures you've posted
    map = [1 1 1;0 0 1;0 1 1;1 1 0;1 0 0];
    
    % Boost histogram values greater than zero so they don't fall in the
    % white band of the colormap.
    S    = size(map,1);
    Hmax = max(max(H));
    bias = (Hmax-S)/(S-1);
    idx = find(H>0);
    H(idx) = H(idx) + bias;
    
    % Plot the histogram
    pcolor([0:Nx-1]*dX+ValMinX, [0:Ny-1]*dY+ValMinY, H);
    shading flat;
    colormap(map);
    

    After increasing N to 1000000, this gives the following plot for the data generated according to your sample:

    16-QAM with AWGN noise - histogram plot