How can I project a MOD04L2 data correctly

I have a MOD04L2 HDF data. I read it and fill NaN in it successfully. When I want to construct coordinates and project the HDF data in it. I find the shape of pixels pulled out of shape. How can I fill the pixels into grids of my coordinates? Can I use sensorA and sensorZ table in the HDF to project my data? My data is in the attachment.

Here is my program

% read data
path = 'MOD04_L2.A2017286.0315.006.2017286133610.hdf'; % my path 
swath_name='mod04';                                    % swath name
lonfield='Longitude';                                  % longitude field name
latfield='Latitude';                                   % latitude field name
datafield='AOD_550_Dark_Target_Deep_Blue_Combined';    % data table field name
flagfield='AOD_550_Dark_Target_Deep_Blue_Combined_QA_Flag';   % flag table field name
lon = hdfread(path,swath_name,'field',lonfield);       % read longitude
lat = hdfread(path,swath_name,'field',latfield);       % read latitude
data = hdfread(path,swath_name,'field',datafield);     % read data table
flag = hdfread(path,swath_name,'field',flagfield);     % read flag table,flag table show the quality of cells
% fill NaN and amend data
data = double(data);               % turn data into double to fill NaN
isFillValue = find(data==-9999);   % find -9999
data(isFillValue)=nan;             % fill NaN in -999
isGood = find(flag<=1);            % find cell whose flag is 0 or 1
data(isGood) = nan;                % fill NaN in low quality cells
data = 0.0010000000474974513*data; % multiple scale factor to calibrate my data
% show my image, I expand my data by 0.5 degrees to show the whole zone
minlon = roundn(min(lon,[],'all'),-1)-0.5;     % min longitude
maxlon = roundn(max(lon,[],'all'),-1)+0.5;     % max longitude
minlat = roundn(min(lat,[],'all'),-1)-0.5;     % min latitude
maxlat = roundn(max(lat,[],'all'),-1)+0.5;     % max latitude
lonlim = [minlon maxlon];                      % longitude range
latlim = [minlat maxlat];                      % latitude range

rastersize = size(data);  %[(maxlat-minlat)*10 (maxlon-minlon)*10];   % I use the size of table as the rastersize, is this phrase wrong?
R = georefcells(latlim,lonlim,rastersize);              % build up R 
[X,Y] = geographicToIntrinsic(R,lat,lon);               % coordinates of cells in my new grid
% demonstrate my image
pcolor(X, Y, data); 
shading faceted; 
colormap("jet");
grid minor
clim([0.0 1.8])
c = colorbar('eastoutside');
c.Label.String = 'AOD';

[Ganesh]

Hi @晨,

What do you mean by "the pixels should fit the grids of my coordinates"?

Do you expect to see a flat image where the coordinates of LAT and LON are directly translated to X and Y values?