% This script reads in data from files generated by NGA
function data = NGAdatareader_large(filename,varidx)

% clear
% clc
% close all
PRODRATES_FLAG = 0;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% User modified variables

%filename = '~/Research/H2_Example/1D_Flame_Detailed/phi0_4/thermal_diffusion_comparison/thermal_diffusion/optdata_1.000E-01'

Lx      = 0.048;
Ly      = 1;
Lz      = 1;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

fid     = fopen(filename,'r');
dims    = fread(fid,4,'integer*4','ieee-le'); % read only the first 4 elements from the file

% Dimensions
nx      = dims(1);
ny      = dims(2);
nz      = dims(3);
nvar    = dims(4);

nsize   = nx*ny*nz;

deltax= Lx/nx;      % mesh size

% The following lines are necessary to move the file pointer to the correct
% location befory you start reading in the variables
fread(fid,1,'real*8','ieee-le'); % skipping over 'dt'
fread(fid,1,'real*8','ieee-le'); % skipping over 'time'
varnames = [];
for var=1:nvar
    varnames = [varnames; fread(fid,8,'*char','ieee-le')]; % skipping over the stored names
end

varnames;

% There are usually 5 variables (or nvar number of variables) stored in the data
% file: U, V, W, P, ZMIX. Each fread call makes the file pointer shift to the end
% of the number of values read. Hence each subsequent fread call will start after
% where the previous fread call stopped. Add on extra variables if necessary.
% %

% if isempty(strfind(filename,'opt'))
%     dummy   = fread(fid,nsize,'real*8','ieee-le');    % will produce a column vector
%     U       = reshape(dummy,nx,ny,nz);                % now turning the column vector into a 3D matrix
%     dummy   = fread(fid,nsize,'real*8','ieee-le');
%     V       = reshape(dummy,nx,ny,nz);
%     dummy   = fread(fid,nsize,'real*8','ieee-le');
%     W       = reshape(dummy,nx,ny,nz);
%     dummy   = fread(fid,nsize,'real*8','ieee-le');
%     P       = reshape(dummy,nx,ny,nz);
%     dummy   = fread(fid,nsize,'real*8','ieee-le');
%     Y{1}    = reshape(dummy,nx,ny,nz);
%     dummy   = fread(fid,nsize,'real*8','ieee-le');
%     Y{2}    = reshape(dummy,nx,ny,nz);
%     
%     data = Y;
% else
%     dummy   = fread(fid,nsize,'real*8','ieee-le');    % will produce a column vector
%     SRC_PROG       = reshape(dummy,nx,ny,nz);                % SRC_PROG
%     
%     data = SRC_PROG;
% end
for i=1:nvar
    dummy   = fread(fid,nsize,'real*8','ieee-le');    % will produce a column vector
    if i == varidx
        data = reshape(dummy,nx,ny,nz);
    end
end
fclose(fid);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Now use U, V, W, P, ZMIX for computations you need
%end