Ocean Color Forum - Not logged in
Forum Ocean Color Home Help Search Login
Previous Next Up Topic SeaDAS / Non-SeaDAS Packages (e.g. MATLAB, ENVI, GIS, etc) / Seawifs to Surfer and Arcview (locked) (6570 hits)
By osinaldi Date 2005-07-14 20:32
I wrote a routine for read the hdf from matlab and write in ARC ASCII grid format and Surfer Grid, basically this routine interpolates the values reported from the variables latitude,longitude and the parameter(my case Chla and k490) and generates a interpolated image in a desidered resolution (my case 0.02 degrees), I want to know if It's valid appears to be, I checked the values with the image in seadas and ok, good but want a opinion.
By osinaldi Date 2005-07-14 20:34
Here's the code, note that it's my draft(a madness) I need to clean, simplificate  and coment It.I for MATLAB v7 (sinceriliy don't know if work in old versions like 6). Converts all the hdf L2 images in the current working directory.

%function [Status] = seawifs2arc
filelist=ls('*.hdf');
filelist=cellstr(filelist);
[h j]=size(filelist);
fileoutc=strrep(filelist,'.hdf','_c.asc');
fileoutk=strrep(filelist,'.hdf','_k.asc');
    to=0.02;
    tu=-99999;
for u =1:h
    chlor_a = hdfread(char(filelist(u,1)),'chlor_a');
    latitude = hdfread(char(filelist(u,1)),'latitude');
    longitude = hdfread(char(filelist(u,1)),'longitude');
    k_490 = hdfread(char(filelist(u,1)),'K_490');
    size(chlor_a);
    [y x]=size(chlor_a);
    nelements=x*y;
    ptsfile=zeros(nelements,4);
    ptsfile(:,1)=reshape(longitude',nelements,1);
    ptsfile(:,2)=reshape(latitude',nelements,1);
    ptsfile(:,3)=reshape(chlor_a',nelements,1);
    ptsfile(:,4)=reshape(k_490',nelements,1);
    latitudes=sort(ptsfile(:,2));
    longitudes=sort(ptsfile(:,1));
    minlat=latitudes(1,1);
    maxlat=latitudes((x*y),1);
    minlong=longitudes(1,1);
    maxlong=longitudes((x*y),1);
    dimx=maxlong-minlong;
    dimy=maxlat-minlat;
    xpts=minlong:.02:maxlong;
    ypts=maxlat:-.02:minlat;
    [Xpts,Ypts] = meshgrid(xpts,ypts);
    Zchlora = griddata(ptsfile(:,1),ptsfile(:,2),ptsfile(:,3),Xpts,Ypts,'linear');
    ZK490 = griddata(ptsfile(:,1),ptsfile(:,2),ptsfile(:,4),Xpts,Ypts,'linear');
    [g h]=size(Zchlora);
    Zchloraz=Zchlora;
    ZK490z=ZK490;
    Zchloraz(isnan(Zchloraz)) = 0;
    ZK490z(isnan(ZK490z)) = 0;
    clorv=sort(reshape(Zchloraz',g*h,1));
    k490v=sort(reshape(ZK490z',g*h,1));
    headerc=zeros(5,2);
    headerc=double(headerc);
    headerc(2,1)=h;
    headerc(2,2)=g;
    headerc(3,1)=minlong;
    headerc(3,2)=maxlong;
    headerc(4,1)=minlat;
    headerc(4,2)=maxlat;
    headerc(5,1)=clorv(1,1);
    headerc(5,2)=clorv(g*h,1);
    headerk=headerc;
    headerk(5,1)=k490v(1,1);
    headerk(5,2)=k490v(g*h,1);
    Zchlora(isnan(Zchlora)) = -99999;
    ZK490(isnan(ZK490)) = -99999;
    %fid=fopen(char(fileoutc(u,1)),'w');
    %fprintf(fid,'DSAA\n',headerc);
    %fclose(fid);
    %dlmwrite(char(fileoutc(u,1)),headerc(2,:-),'delimiter',' ','precision','%0.0f','-append');
    fid=fopen(char(fileoutc(u,1)),'w');
    fprintf(fid,'NCOLS %g \n',headerc(2,1));
    fprintf(fid,'NROWS %g \n',headerc(2,2));
    fprintf(fid,'XLLCORNER %g \n',headerc(3,1));
    fprintf(fid,'YLLCORNER %g \n',headerc(4,1));
    fprintf(fid,'CELLSIZE %g \n',to);
    fprintf(fid,'NODATA_VALUE %g \n',tu);
    fclose(fid)
    dlmwrite(char(fileoutc(u,1)),Zchlora,'delimiter',' ','precision','%10.10f','-append');
    fid=fopen(char(fileoutk(u,1)),'w');
    fprintf(fid,'NCOLS %g \n',headerc(2,1));
    fprintf(fid,'NROWS %g \n',headerc(2,2));
    fprintf(fid,'XLLCORNER %g \n',headerc(3,1));
    fprintf(fid,'YLLCORNER %g \n',headerc(4,1));
    fprintf(fid,'CELLSIZE %g \n',to);
    fprintf(fid,'NODATA_VALUE %g \n',tu);
    fclose(fid)
    dlmwrite(char(fileoutk(u,1)),Zchlora,'delimiter',' ','precision','%10.10f','-append');
    tprocd=clock;
    fprintf('Image done: %s number %g of %g at %g:%g:%g\n',char(filelist(1,1)),u,h,tprocd(1,4),tprocd(1,5),tprocd(1,6));
    Status='Done!'
    end
   
By Anonymous Date 2006-05-23 22:46
I need a routine for read chla -a oc4 seawifs hdf from matlab, basically interpolates tha values reported from latitude and longitude and the parameter (chla-a oc4) and generates a interpolated image in a desired resolution. I new in matlab please help me to write m-file (i am using m-file of osinaldi but not work (Error in ==> hdfread at 219; [hinfo,subsets] = dataSetInfo(varargin{:});Error in ==> seawifs2arc at 10; chlor_a = hdfread(char(filelist(u,1)),'chlor_a');-).

The m file complete is:

function [Status] = seawifs2arc
filelist=ls('*.hdf');
filelist=cellstr(filelist);
[h j]=size(filelist);
fileoutc=strrep(filelist,'.hdf','_c.asc');
fileoutk=strrep(filelist,'.hdf','_k.asc');
    to=0.02;
    tu=-99999;
for u =1:h
    chlor_a = hdfread(char(filelist(u,1)),'chlor_a');
    latitude = hdfread(char(filelist(u,1)),'latitude');
    longitude = hdfread(char(filelist(u,1)),'longitude');
    size(chlor_a);
    [y x]=size(chlor_a);
    nelements=x*y;
    ptsfile=zeros(nelements,4);
    ptsfile(:,1)=reshape(longitude',nelements,1);
    ptsfile(:,2)=reshape(latitude',nelements,1);
    ptsfile(:,3)=reshape(chlor_a',nelements,1);
    latitudes=sort(ptsfile(:,2));
    longitudes=sort(ptsfile(:,1));
    minlat=latitudes(1,1);
    maxlat=latitudes((x*y),1);
    minlong=longitudes(1,1);
    maxlong=longitudes((x*y),1);
    dimx=maxlong-minlong;
    dimy=maxlat-minlat;
    xpts=minlong:.02:maxlong;
    ypts=maxlat:-.02:minlat;
    [Xpts,Ypts] = meshgrid(xpts,ypts);
    Zchlora = griddata(ptsfile(:,1),ptsfile(:,2),ptsfile(:,3),Xpts,Ypts,'linear');
    ZK490 = griddata(ptsfile(:,1),ptsfile(:,2),ptsfile(:,4),Xpts,Ypts,'linear');
    [g h]=size(Zchlora);
    Zchloraz=Zchlora;
    Zchloraz(isnan(Zchloraz)) = 0;
    clorv=sort(reshape(Zchloraz',g*h,1));
    headerc=zeros(5,2);
    headerc=double(headerc);
    headerc(2,1)=h;
    headerc(2,2)=g;
    headerc(3,1)=minlong;
    headerc(3,2)=maxlong;
    headerc(4,1)=minlat;
    headerc(4,2)=maxlat;
    headerc(5,1)=clorv(1,1);
    headerc(5,2)=clorv(g*h,1);
    headerk=headerc;
    Zchlora(isnan(Zchlora)) = -99999;
    %fid=fopen(char(fileoutc(u,1)),'w');
    %fprintf(fid,'DSAA\n',headerc);
    %fclose(fid);
    %dlmwrite(char(fileoutc(u,1)),headerc(2,:-),'delimiter',' ','precision','%0.0f','-append');
    fid=fopen(char(fileoutc(u,1)),'w');
    fprintf(fid,'NCOLS %g \n',headerc(2,1));
    fprintf(fid,'NROWS %g \n',headerc(2,2));
    fprintf(fid,'XLLCORNER %g \n',headerc(3,1));
    fprintf(fid,'YLLCORNER %g \n',headerc(4,1));
    fprintf(fid,'CELLSIZE %g \n',to);
    fprintf(fid,'NODATA_VALUE %g \n',tu);
    fclose(fid)
    dlmwrite(char(fileoutc(u,1)),Zchlora,'delimiter',' ','precision','%10.10f','-append');
    fid=fopen(char(fileoutk(u,1)),'w');
    fprintf(fid,'NCOLS %g \n',headerc(2,1));
    fprintf(fid,'NROWS %g \n',headerc(2,2));
    fprintf(fid,'XLLCORNER %g \n',headerc(3,1));
    fprintf(fid,'YLLCORNER %g \n',headerc(4,1));
    fprintf(fid,'CELLSIZE %g \n',to);
    fprintf(fid,'NODATA_VALUE %g \n',tu);
    fclose(fid)
    dlmwrite(char(fileoutk(u,1)),Zchlora,'delimiter',' ','precision','%10.10f','-append');
    tprocd=clock;
    fprintf('Image done: %s number %g of %g at %g:%g:%g\n',char(filelist(1,1)),u,h,tprocd(1,4),tprocd(1,5),tprocd(1,6));
    Status='Done!'
    end
Previous Next Up Topic SeaDAS / Non-SeaDAS Packages (e.g. MATLAB, ENVI, GIS, etc) / Seawifs to Surfer and Arcview (locked) (6570 hits)



Responsible NASA Official: Gene C. Feldman
Curator: OceanColor Webmaster
Authorized by: Gene C. Feldman
Updated: 27 November 2007 		Privacy Policy and Important Notices NASA logo