Extract data from specific latitudes and longitudes in time in NetCDF (.nc) to .csv

I am using Matlab 2021a on Linux Mint. I am trying to extract data from a file in NetCDF, from 3 specific latitude and 3 specific longitudes for the whole time series. At the time of generating a .csv file, it is giving an error in the 'writematrix' command. I tried another way with the 'xlswrite' command and was not successful.

The link to download the CPC rainfall data I inserted just below in the script.

I am very grateful for the answer.

Assinado,

Augusto Pereira.

% Download the CPC data used in the script below
% https://psl.noaa.gov/thredds/catalog/Datasets/cpc_global_precip/catalog.html
% precip.2021.nc


LG=[-48.25 -48.75 -49.25];
LT=[-1.25 -1.25 -1.25];

filename='cpc_global_precip_precip.2021.nc'
ncdisp(filename,'/','min');
%ncdisp(filename)

%precip='pr'
Precip=ncread(filename,'precip');
long=ncread(filename,'lon');
lat=ncread(filename,'lat');
time=ncread(filename,'time');

for i=1:3 
LGG=find(long==LG(i));
LTT=find(lat==LT(i));
 if isempty(LGG)||isempty(LTT)
    disp('LOCATION NOT FOUND IN THIS NetCDF FILE')
    break
 end
% Step 5 :Convert 3 dimensional array to column vector array for one station
Precip(:,i)=precip(LGG,LTT,:);
X(i)={'Precipitation(mm/day)'};% label for variable
end
latt2=[LG ; LT];

time=double(time);
AA=time/24+datenum('1900-01-01 00:00:0.0');
yy=year(AA);mm=month(AA);dd=day(AA);
date1=[yy mm dd];
DD={'Year','Month','Day'};

filename1='PRP_CPC.csv';

sheet=1;
C={'LOCATION: Cities','','','Longitude';'DATE FROM JAN-DEZ 2021','','','Latitude'};


writematrix(filename1,C,sheet,'A1');
writematrix(filename1,DD,sheet,'A3');
writematrix(filename1,latt2,sheet,'E1');
writematrix(filename1,X,sheet,'E3');
writematrix(filename1,date1,sheet,'A4');
writematrix(filename1,Precip,sheet,'E4');

[Nicholas Darangwa]

I am also having this same challenge, extracting climate data from an netcdf file, need time series for given set of coodinates from file on this link:

I had tried work on the above code to this point, but getting errors.

clear all,clc

% Step 1: Specify longtitude and and latitude

LG=[16.5 19.5];

LT=[-33 -35];

% Step 2: Display all variables in NetCDF file using dcsisplay command

filename = 'cams73_latest_ch4_conc_surface_satellite_dm_200903.nc';

ncdisp(filename,'/','min');

% Step 3: Read variables from NetCDF file by command ncread

Meth_Levels = ncread(filename,'CH4');% Dimensions: lon,lat,time

long = ncread(filename,'longitude');

lat = ncread(filename,'latitude');

time = ncread(filename,'time');

for i=1:3

LGG=find(long==LG(i));

LTT=find(lat==LT(i));

if isempty(LGG)||isempty(LTT)

disp('LOCATION NOT FOUND IN THIS NetCDF FILE')

break

end

% Step 5 :Convert 3 dimensional array to column vector array for one station

Meth_Levels(:,i)=Meth_Levels(LGG,LTT,:); % Unrecognized function or variable 'precip'.

X(i)={'Methane_Data (mm/day)'};% label for variable

end

latt2=[LG ; LT];

time=double(time);

AA=time/24+datenum('1900-01-01 00:00:0.0');

[yy,mm,dd,~,~,~] = datevec(AA);

date1=[yy mm dd];

DD={'Year','Month','Day'};

The error says:

Unable to perform assignment because the size of the left side is 120-by-1 and the size of the right side is

1-by-1-by-1054.

[Nicholas Darangwa]

Here is the link for the file: NetCDF

[Nicholas Darangwa]

The intention is to extract the data from multiple files, was working on this, and testing, so i have this error that i dont know how to fix...

[Mathieu NOE]

hello

in the future please create a new post so that we can better follow / answer

anyway , here is the solution

NB that the CH4 data is 4D and not 3D,

also do not use the same variable name on the LHS and RHS of an equation, especially if the dimensions are not the same on both side

the for loop cannot run 3 iterations as your input vectors LG / LT has only 2 values (minor bug)

try this now :

clear all,clc
 
% Step 1: Specify longtitude and and latitude
 
LG=[16.5 19.5];
LT=[-33 -35];
 
% Step 2: Display all variables in NetCDF file using dcsisplay command
 
filename = 'cams73_latest_ch4_conc_surface_satellite_dm_200903.nc';
ncdisp(filename,'/','min');
 
% Step 3: Read variables from NetCDF file by command ncread
 
Meth_Levels = ncread(filename,'CH4');
%            Size:       120x90x34x31
%            Dimensions: longitude,latitude,level,time

long = ncread(filename,'longitude');
lat = ncread(filename,'latitude');
time = ncread(filename,'time');

% return
 
for i=1:numel(LG) 
    LGG=find(long==LG(i));
    LTT=find(lat==LT(i));
     if isempty(LGG)||isempty(LTT)
        disp('LOCATION NOT FOUND IN THIS NetCDF FILE')
        break
     end
    % Step 5 :Convert 4 dimensional array to column vector array for one station
    M_Levels(:,:,i) = squeeze(Meth_Levels(LGG,LTT,:,:));  % do not use the same variable name on the LHS and RHS
    X(i)={'Methane_Data (mm/day)'};% label for variable
end
latt2=[LG ; LT];
time=double(time);
AA=time/24+datenum('1900-01-01 00:00:0.0');
[yy,mm,dd,~,~,~] = datevec(AA);
date1=[yy mm dd];
DD={'Year','Month','Day'};

[Nicholas Darangwa]

Thank you so much. So much appreciated. This one runs with no errors.

[Mathieu NOE]

my pleasure !

Answer 1

[Mathieu NOE]

hello Augusto

I fixed a few minor bugs in the code .

Also I noticed from the nc file that the range for longitude data is 0 to 360° and not -180 / + 180° so if you're asking for a negative long point it will return empty answer. I had to add 360 ° to your negative input values to make it coherent with the file

The nc file I downloaded has a different name but I beleive it's the same as yours

see attached output file for info

Updated code :

% Download the CPC data used in the script below
% https://psl.noaa.gov/thredds/catalog/Datasets/cpc_global_precip/catalog.html
% precip.2021.nc
LG=[-48.25 -48.75 -49.25] + 360;  % longitude (NB range is 0 : 360°)
LT=[-1.25 -1.25 -1.25];     % latitude
% filename='cpc_global_precip_precip.2021.nc'
filename='precip.2021.nc'
ncdisp(filename,'/','min');
precip=ncread(filename,'precip');  % Dimensions: lon,lat,time
long=ncread(filename,'lon');
lat=ncread(filename,'lat');
time=ncread(filename,'time');
for i=1:3 
    LGG=find(long==LG(i));
    LTT=find(lat==LT(i));
     if isempty(LGG)||isempty(LTT)
        disp('LOCATION NOT FOUND IN THIS NetCDF FILE')
        break
     end
    % Step 5 :Convert 3 dimensional array to column vector array for one station
    Precip(:,i)=precip(LGG,LTT,:);  % Unrecognized function or variable 'precip'.
    X(i)={'Precipitation(mm/day)'};% label for variable
end
latt2=[LG ; LT];
time=double(time);
AA=time/24+datenum('1900-01-01 00:00:0.0');
[yy,mm,dd,~,~,~] = datevec(AA);
date1=[yy mm dd];
DD={'Year','Month','Day'};
filename1='PRP_CPC.xlsx';
sheet=1;
C={'LOCATION: Cities','','','Longitude';'DATE FROM JAN-DEZ 2021','','','Latitude'};
writecell(C,filename1,"Sheet",sheet,"Range",'A1');
writecell(DD,filename1,"Sheet",sheet,"Range",'A3');
writematrix(latt2,filename1,"Sheet",sheet,"Range",'E1');
writecell(X,filename1,"Sheet",sheet,"Range",'E3');
writematrix(date1,filename1,"Sheet",sheet,"Range",'A4');
writematrix(Precip,filename1,"Sheet",sheet,"Range",'E4');

Answer 2

[Augusto Gabriel da Costa Pereira]

Thank you very much for the answer, Mathieu. Instead of saving in .xlsx via writematrix, can I save in .mat matrix? If so, how? I am using linux, so there is extension error when saving in .xlsx.

[Mathieu NOE]

hello Augusto

have you asked TheMathWorks for support regarding this issue ?

yes you can save one or several variablesto mat file but how you want then to use it afterwards ? I suppose you want anyway an excel file at the end

Save workspace variables to file - MATLAB save - MathWorks France

[Augusto Gabriel da Costa Pereira]

My teachers only work with matlab, so I wanted to save in .mat, the Excel file would stay with me just for security purposes.

[Mathieu NOE]

ok

did you read the doc for save ?

Save workspace variables to file - MATLAB save - MathWorks France

[Augusto Gabriel da Costa Pereira]

I managed to save. Thanks a lot for the help.

[Augusto Gabriel da Costa Pereira]

Hello Matheus. I have another doubt.

That script you fixed for me, I can extract data from just one NetCDF file.

If I want to extract data from multiple NetCDF for a specific latitude and longitude, how can I do that? If yes, do you have the script for that?

[Mathieu NOE]

hello again

the code below will list all .nc files in the selected directory and do the data extraction and save to output file in a (sub)function. I kept the excel output but you can change that to mat as you already did;

hope it helps

LG=[-48.25 -48.75 -49.25] + 360;  % longitude (NB range is 0 : 360°)
LT=[-1.25 -1.25 -1.25];     % latitude

%% main code 

fileDir = pwd; % current directory (or specify which one is the working directory)

S = dir(fullfile(fileDir,'*.nc')); % get list of nc files in directory
S = natsortfiles(S); % sort file names into natural order , see : 
           %(https://fr.mathworks.com/matlabcentral/fileexchange/47434-natural-order-filename-sort)
for k = 1:length(S)
    filename = S(k).name % to actually show filenames are sorted (see command window)
    sheet = 1; % for excel file output, either on first sheet or increment sheet number if needed
    filename_out = [filename(1:length(filename)-3) '.xlsx']; % replace nc extension with xlsx (or whatever you want)
    get_and_save_ncf_data(fullfile(fileDir, filename),LG,LT,fullfile(fileDir, filename_out),sheet);

end

%%%%%%%%%%%%%%%% functions %%%%%%%%%%%%%%%%%%%%%%%

function get_and_save_ncf_data(filename,LG,LT,filename_out,sheet)

% do all the recurring tasks in a separate function below 

%     ncdisp(filename,'/','min');
    precip=ncread(filename,'precip');  % Dimensions: lon,lat,time
    long=ncread(filename,'lon');
    lat=ncread(filename,'lat');
    time=ncread(filename,'time');
    for i=1:3 
        LGG=find(long==LG(i));
        LTT=find(lat==LT(i));
         if isempty(LGG)||isempty(LTT)
            disp('LOCATION NOT FOUND IN THIS NetCDF FILE')
            break
         end
        % Step 5 :Convert 3 dimensional array to column vector array for one station
        Precip(:,i)=precip(LGG,LTT,:);  % Unrecognized function or variable 'precip'.
        X(i)={'Precipitation(mm/day)'};% label for variable
    end
    latt2=[LG ; LT];
    time=double(time);
    AA=time/24+datenum('1900-01-01 00:00:0.0');
    [yy,mm,dd,~,~,~] = datevec(AA);
    date1=[yy mm dd];
    DD={'Year','Month','Day'};
    C={'LOCATION: Cities','','','Longitude';'DATE FROM JAN-DEZ 2021','','','Latitude'};
    % export to excel
    writecell(C,filename_out,"Sheet",sheet,"Range",'A1');
    writecell(DD,filename_out,"Sheet",sheet,"Range",'A3');
    writematrix(latt2,filename_out,"Sheet",sheet,"Range",'E1');
    writecell(X,filename_out,"Sheet",sheet,"Range",'E3');
    writematrix(date1,filename_out,"Sheet",sheet,"Range",'A4');
    writematrix(Precip,filename_out,"Sheet",sheet,"Range",'E4');
    disp(['Data exported to file :' filename_out])
end