# How can I georeference this grid/matrix?

3 Ansichten (letzte 30 Tage)
Harith Kubaisy am 26 Jul. 2020
Kommentiert: Maskus Kit am 6 Feb. 2021
Hello,
I have a script which plots kernel density values for earthquakes. However, I need to visualize the final results of the p array on a proper projected geograhpic map. How can I incorporate the results on a georeferenced MATLAB map? If that is not possible, what is the best way to export the data to ArcGIS platform?
I tried saving the results as a .txt file (ASCII) and performed a conversion to Raster by ArcMap 10.4.1 but I could not get meaningful results.
The script can be seen below, and the longitude/latitude files are attached.
clc,clear
long = (long5);
lat = (lat5);
step = 0.009;
x = min(long):step:max(long);
y = min(lat):step:max(lat);
[X,Y]= ndgrid(x,y); % creates grid
R = 0.05; % circle radius
plot(long,lat,'.r') % plot data
P = X*0; % no. of points inside each circle
T = 2019-1916+1;
sig = 2;
hold on
for i= 1:length(x)
for j=1:length(y)
D2 = (long-x(i)).^2 + (lat-y(j)).^2;
ix = D2 < R^2; % find smaller distances
D(ix) = scale*D2(ix); % write smaller distances
P(i,j) = sum(exp(-D(ix)/2/sig^2))/T; %write ID
N(i,j) = sum(ix); % no. of points inside (i,j) circle
% plot(long(ix),lat(ix),'ob') %plot data inside circle (if exists)
if sum(ix) % if there are points inside
viscircles([x(i),y(j)], R,'edgeColor','g'); % plot circle
else
% viscircles([x(i),y(j)], R,'edgecolor','r'); % plot circle
end
end
end
p = P/(pi*6.^2);
G = rot90(flip(p),3);
imagesc(x,y,G);
colorbar
hold off
axis equal
##### 3 Kommentare1 älteren Kommentar anzeigen1 älteren Kommentar ausblenden
Harith Kubaisy am 5 Feb. 2021
Hello. Yes the following work-around did it for me. Sorry for not updating the post.
In a seprate script:
V1 = X(:);
V2 = Y(:);
PV = pMo(:);
lon = (V1);
lat = (V2);
weights = (PV);
hold off
figure
geodensityplot(lat,lon,weights,'FaceColor','interp');
hold on
geobasemap grayterrain
colormap(hsv)
colorbar
Maskus Kit am 6 Feb. 2021

Melden Sie sich an, um zu kommentieren.

### Antworten (1)

jonas am 26 Jul. 2020
With the mapping toolbox you can just replace your plot functions with the corresponding map axes function.
ax = worldmap('italy');
plotm(coastlat,coastlon)
plotm(long,lat,'.r') % plot data
h = pcolorm(x,y,G);
should give you something like this
##### 1 Kommentar-1 ältere Kommentare anzeigen-1 ältere Kommentare ausblenden
Harith Kubaisy am 27 Jul. 2020
Thank you for the quick response, I tried the approach but I did not get the plot over my desired area.
See attached script and result map.

Melden Sie sich an, um zu kommentieren.

### Kategorien

Mehr zu Mapping Toolbox finden Sie in Help Center und File Exchange

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by