data files processing and plotting
1 view (last 30 days)
Show older comments
I have two years dataset which is in daily format. I have following codes to read this data
path=('/datadir/gsat/');
co2_month=[];
lat_month=[];
lon_month=[];
files=dir('*h5');
monthsStr={'01';'02';'03';'04';'05';'06';'07';'08';'09';'10';'11';'12'};
monstring={'January','February','March','April','May','June','July','August','September','October','November','December'};
for ii=1:numel(files)
for year=2009:2010
for jj=1:numel(monthsStr)
nameFound=strfind(files(ii).name,['GOSATTFTS' num2str(year) monthsStr{jj}]);
if nameFound
latitude=hdf5read(strcat(path,files(ii).name),'/Data/geolocation/latitude');
longitude=double(hdf5read(strcat(path,files(ii).name),'/Data/geolocation/longitude'));
xco2=double(hdf5read(strcat(path,files(ii).name),'/Data/mixingRatio/XCO2'));
time =hdf5read(strcat(path,files(ii).name),'/scanAttribute/time');
times=length(latitude);
% find lat, long for area of interest
d_deg=1.5;
match=find(latitude < (45+d_deg/2) & latitude >=(15-d_deg/2) & longitude < (-75)+d_deg/2 & longitude >=(-105)-d_deg/2);
lat_G=double(latitude(match));
lon_G=double(longitude(match));
xco2_G=double(xco2(match));
if ~isempty(xco2_G)
co2_month=[co2_month;xco2_G];
lat_month=[lat_month;lat_G];
lon_month=[lon_month;lon_G];
end
d_deg=0.5;
latgrid=15:d_deg:45;
longrid=(-105):d_deg:(-75);
latind=round(lat_month/d_deg-15/d_deg+1);
lonind=round(lon_month/d_deg -(-105/d_deg)+1);
sum_co2=zeros(length(latgrid),length(longrid));
cnt_co2=zeros(length(latgrid),length(longrid));
for j=1:length(co2_month)
sum_co2(latind(j),lonind(j))=sum_co2(latind(j),lonind(j))+co2_month(j);
cnt_co2(latind(j),lonind(j))=cnt_co2(latind(j),lonind(j))+1;
end
grid_co2=sum_co2./cnt_co2;
end
end
end
end
when ever I compile this codes it reads all data .
But I would like to plot every months data on one map e.g observations from 1-31 days on January map. I am using following plotting codes
figure(1)
load coast
h1=axesm('MapProjection','eqdcylin',...
'Grid','off',...
'MapLatLimit',[15-d_deg/2 45+d_deg/2],...
'MapLonLimit',[(-105-d_deg/2) (-75+d_deg/2)],...
'flatlimit',[15-d_deg/2 45+d_deg/2],...
'flonlimit',[(-105-d_deg/2) (-75+d_deg/2)],...
'MeridianLabel','off',...
'Frame','on');
plotm(lat,long,'k')
hold on
tightmap on
colormap(Jet);
grid on
h=surfacem(latgrid,longrid,grid_co2);
so somebody please tell me how I can do it in matlab
0 Comments
Answers (1)
Fangjun Jiang
on 19 Aug 2011
You need to think about the flow of your code. If you want to plot the data of every month into one plot, you need to go through a for-loop iterating on the month.
Maybe start with a simple one, pick a month, find all the data belongs to that month, load the data and plot them into one plot. Then, wrap that code inside a for-loop.
Your current code is not efficient. For every particular year and month, you have to compare it to all those more than 700 file names every time. That is waste of time. Once you determine the year and month, it's easy to pick out the corresponding 30 files. Please go back to look at my answer to your previous question. It should not be that hard to wrap that solution into a for-loop to achieve what you try to do here.
0 Comments
See Also
Categories
Find more on Map Display in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!