Showing posts with label offset. Show all posts
Showing posts with label offset. Show all posts

Saturday, February 25, 2012

[Matlab] show the variable information from a netcdf file

function myShowNc(ncFileName)
% show information of a given netcdf file, similar to ncdump
%          using matlab native netcdf functions
% usage:
%           myShowNc(netcdf-filename)
% http://scriptdemo.blogspot.com

if nargin==0
   help myShowNc
   return
end
if ~exist(ncFileName,'file')
   error([ncFileName,' is not found!']);
end

ncfid=netcdf.open(ncFileName,'NC_NOWRITE');

% check dimensions
[ndims,nvars,ngatts,unlimdimID] = netcdf.inq(ncfid);
myncinfo.unLimID=unlimdimID;
myncinfo.dimname=cell(ndims,1);
myncinfo.dimlength=zeros(ndims,1);
myncinfo.varname=cell(nvars,1);
myncinfo.vardims=cell(nvars,1);
myncinfo.varscale=ones(nvars,1);
myncinfo.varoffset=zeros(nvars,1);

disp(['netcdf file: ',ncFileName])
disp(' dimension:')
for Ndim=1:ndims
    [myncinfo.dimname{Ndim},myncinfo.dimlength(Ndim)]=netcdf.inqDim(ncfid,Ndim-1);
    disp([' ',fixStrLen(myncinfo.dimname{Ndim},10),': ',num2str(myncinfo.dimlength(Ndim))])
end

disp(' variables:')
disp([' ',fixStrLen('var-name |',12,0),fixStrLen('var-dims',30,0),fixStrLen('scale',10,1),fixStrLen('offset',10,1)])
% check variables
for nv=1:nvars
    [myncinfo.varname{nv},myncinfo.vardims{nv}, ...
     myncinfo.varscale(nv),myncinfo.varoffset(nv)]=getVarInfo(ncfid,nv-1,myncinfo);
    disp([' ', fixStrLen([myncinfo.varname{nv},' |'],12) , fixStrLen(myncinfo.vardims{nv},30), ...
                fixStrLen(num2str(myncinfo.varscale(nv)),10,1),fixStrLen(num2str(myncinfo.varoffset(nv)),10,1)])
end
netcdf.close(ncfid);

%--------------------------------------------------------------------------------------
function outStr=fixStrLen(inStr,nLen,opt)
% renamed to formatStr in another post.
% if  opt == 0  'right-align' [default case]
%     opt == 1  'center-align'
%     opt == 2  'left-align'
if nargin==1
   nLen=0;
elseif nargin==2
   opt=0;
elseif nargin~=3
   help fixStrLen;
   return
end

oriLen=length(inStr);
if (oriLen>=nLen)
    outStr=inStr;
else
    switch (opt)
        case {0}
            nSpace=nLen-oriLen;
            outStr=[repmat(' ',1,nSpace),inStr];
        case {1}
          nSpace=ceil(0.5*(nLen-oriLen));
          outStr=[repmat(' ',1,nSpace),inStr];
          nSpace=nLen-oriLen-nSpace;
          if (nSpace>0)
             outStr=[outStr,repmat(' ',1,nSpace)];
          end
        case {2}
            nSpace=nLen-oriLen;
            outStr=[inStr,repmat(' ',1,nSpace)];
        otherwise
            disp('Unknown option');
            return
    end
end

function [varname,dimsStr,myScale,myOffSet]=getVarInfo(ncfid,varID,myncinfo)
[varname,vartype,vardimids,numAtts]=netcdf.inqVar(ncfid,varID);
myScale=1; myOffSet=0;
for numA=1:numAtts
    tmpattname=netcdf.inqAttName(ncfid,varID,numA-1);
    switch lower(tmpattname)
        case {'scale_factor','scalefactor'}
           myScale=netcdf.getAtt(ncfid,varID,tmpattname);
        case {'add_offset','off_set','offset','offset_value'}
           myOffSet=netcdf.getAtt(ncfid,varID,tmpattname);
        otherwise
           % do nothing
     end
end
% get dimStr
if ~isempty(vardimids)
   dimsStr=myncinfo.dimname{vardimids(1)+1};
   for nd=2:length(vardimids)
         dimsStr=[myncinfo.dimname{vardimids(nd)+1},' x ',dimsStr];
   end
else
   dimsStr='unlimited';
end 

Monday, January 3, 2011

[UPDATED] read netcdf file using native matlab functions

check the ncread for new version matlab from Mathworks. 
 
The attributes of scale factor, offset value and missing value are considered in this version.

function ncvardata=GetNcVarMatlab(ncfile,ncvar,RangeStr)
% to read a variable from the netcdf file using native matlab netcdf-functions
% Usage: vardata=GetNcVarMatlab(nc-filename,nc-varname,Dimension-Range)
% e.g.,
%        mydata=GetNcVarMatlab('mask.nc','tmask','1:1,1:46,1:400,1:568');
%                      RangeStr here is similar to the output of ncdump
% note: by default, variable will be converted into the class of double,
%                      and the missing value will be converted into nan.
% copyright @ http://scriptdemo.blogspot.com


IsDebug=0;
IsMissingNaN=1; % convert missing values into NaN
if ~exist(ncfile,'file')
   error([ncfile,' can not been found']);
   return
end

ncfid=netcdf.open(ncfile,'NC_NOWRITE');
% get dimension length
[numdims, numvars, numglobalatts, unlimdimID] = netcdf.inq(ncfid);
for NDim=1:numdims
   [dimname{NDim}, dimlen(NDim)] = netcdf.inqDim(ncfid,NDim-1);
end

%var infor
varid=netcdf.inqVarID(ncfid,ncvar);
[varname,varxtype,vardimids,NumAtts]=netcdf.inqVar(ncfid,varid);
numdims=length(vardimids); % some variables may have less dimensions

%check scale_factor & off_set
IsRescale=0; IsMissing=0;
myScale=1; myOffSet=0.0;
for numA=1:NumAtts
    tmpattname=netcdf.inqAttName(ncfid,varid,numA-1);
    switch lower(tmpattname)
       case {'scale_factor','scalefactor'}
            IsRescale=1;
            myscale=netcdf.getAtt(ncfid,varid,tmpattname);
       case {'add_offset','off_set','offset'}
            IsRescale=1;
            myoffset=netcdf.getAtt(ncfid,varid,tmpattname);
       case {'missing_value','missingvalue','missing'}
            IsMissing=1;
            mymissing=netcdf.getAtt(ncfid,varid,tmpattname);
       otherwise
            % do nothing
    end
end
if IsRescale==1
   if (myscale==1 & myoffset==0)
      IsRescale=0;
   end
else
   if ~exist('myscale','var') myscale=1; end
   if ~exist('myoffset','var') myoffset=0; end
end

if ~exist('RangeStr','var')
   RangeStr=':';
   for nd=2:numdims
       RangeStr=[RangeStr,',:'];
   end
end

IndC=strfind(RangeStr,',');% index of commas
if length(IndC)~=numdims-1
   disp(['Input dimension range is not correct, should be: ',num2str(numdims),'-d']);
   return
end

startstr=''; countstr=''; stridestr='';
for nn=1:length(IndC)+1
    if nn==1
       tmpstr=RangeStr(1:IndC(1)-1);
    elseif nn<numdims
       tmpstr=RangeStr(IndC(nn-1)+1:IndC(nn)-1);
    else
       tmpstr=RangeStr(IndC(nn-1)+1:end);
    end

    [tmpstart,tmpcount,tmpstride]=GetIndFromStr(tmpstr,':',dimlen(vardimids(numdims-nn+1)+1));
    if nn==1
        startstr=[num2str(tmpstart),'],'];
        countstr=[num2str(tmpcount),'],'];
        stridestr=[num2str(tmpstride),']'];
    elseif nn==numdims
        startstr=[',[',num2str(tmpstart),',',startstr];
        countstr=['[',num2str(tmpcount),',',countstr];
        stridestr=['[',num2str(tmpstride),',',stridestr];
    else
        startstr=[num2str(tmpstart),',',startstr];
        countstr=[num2str(tmpcount),',',countstr];
        stridestr=[num2str(tmpstride),',',stridestr];
    end
end

if IsDebug
   disp(['ncvardata=netcdf.getVar(ncfid,varid',startstr,countstr,stridestr,');']);
end
eval(['ncvardata=netcdf.getVar(ncfid,varid',startstr,countstr,stridestr,');']);
netcdf.close(ncfid);
ncvardata=permute(ncvardata,length(size(ncvardata)):-1:1);

if IsRescale==1
   %convert all data type into double
   if ~strcmp(class(ncvardata),'double')
      ncvardata=double(ncvardata);
   end
   if IsMissing
      if IsMissingNaN==1; % convert missing values into NaN
         ncvardata(ncvardata==mymissing)=nan;
         ncvardata=ncvardata*double(myscale)+double(myoffset);
      else
         disp(['missing value for ',ncvar,' is ',num2str(mymissing)]);
         %IndValid=find(ncvardata~=mymissing); length(IndValid)
         %ncvardata(IndValid)=ncvardata(IndValid)*double(myscale)+double(myoffset);
         ncvardata(ncvardata~=mymissing)=ncvardata(ncvardata~=mymissing)*double(myscale)+double(myoffset);
      end
   else
      ncvardata=ncvardata*double(myscale)+double(myoffset);
   end
else
  ncvardata=double(ncvardata);
end

function [numstart,numcount,numstride]=GetIndFromStr(tmpstr,ctag,Ndim)
% Read the numbers from the special string
IsDebug=0;
IndC=strfind(tmpstr,ctag);
if isempty(IndC)
   numstart=str2double(tmpstr)-1;
   numcount=1;
   numstride=1;
else
  switch length(IndC)
      case {1}
          switch IndC
            case {1} % : and :num
                if strcmp(tmpstr,':')
                   numstart=0; numcount=Ndim; numstride=1;
                else
                   if strcmp(tmpstr(2:end),'end')
                    numstart=0; numcount=Ndim; numstride=1;
                   else
                    numstart=0; numcount=str2double(tmpstr(2:end))-1; numstride=1;
                   end
                end
            case {length(tmpstr)}
                if strcmp(tmpstr,':')
                   numstart=0; numcount=Ndim; numstride=1;
                   return
                else
                   numstart=str2double(tmpstr(1:IndC-1))-1;
                   numcount=Ndim-numstart;
                   numstride=1;
                end
            otherwise
               numstart=str2double(tmpstr(1:IndC-1))-1;
               if strcmp(tmpstr(IndC+1:end),'end')
                  numend=Ndim;
               else
                  numend=str2double(tmpstr(IndC+1:end));
               end
               numcount=numend-numstart;
               numstride=1;
          end
      case {2}
          % num0:stride:num1
          if IndC(1)==1
             numstart=0;
          else
             numstart=str2double(tmpstr(1:IndC(1)-1))-1;
          end
          if IndC(2)==IndC(1)+1
             numstride=1;
          else
             numstride=str2double(tmpstr(IndC(1)+1:IndC(2)-1));
          end
          if IndC(2)==length(tmpstr)
                 numend=Ndim-1;
          else
            if strcmp(tmpstr(IndC(2)+1:end),'end')
               numend=Ndim-1;
            else
               numend=str2double(tmpstr(IndC(2)+1:end))-1;
            end
          end
          numcount=floor((numend-numstart)/numstride)+1;
          numend=numstart+(numcount-1)*(numStride);
      otherwise
          error('Not defined range type')
  end
end
if IsDebug
   disp(['tmpstr=',tmpstr])
   disp(['start:',num2str(numstart), ' stride:',num2str(numstride),' count:',num2str(numcount)])
end

ShowCalendar