function [Fz Fr Ft z_values] = hfss_to_opal(fname, Nz, Nr, FREQUENCY)

Tz = load(strcat(fname,"/Ez.fld"));
Tr = load(strcat(fname,"/Ex.fld"));
Tt = load(strcat(fname,"/Hy.fld"));
z_values = unique(Tz(:,3));
r_values = unique(Tz(:,1));
Tz = Tz(:,4:5);
Tr = Tr(:,4:5);
Tt = Tt(:,4:5);
Tz(isnan(Tz(:))) = 0.0;
Tr(isnan(Tr(:))) = 0.0;
Tt(isnan(Tt(:))) = 0.0;

Nz = min(size(z_values,1), Nz);
Nr = min(size(r_values,1), Nr);
Fz = zeros(Nr,Nz);
Fr = zeros(Nr,Nz);
Ft = zeros(Nr,Nz);

for i = 1:Nr
  Fz(i,:) = Tz((i-1) * Nz + 1:i*Nz, 1);
  Fr(i,:) = Tr((i-1) * Nz + 1:i*Nz, 1);
  Ft(i,:) = Tt((i-1) * Nz + 1:i*Nz, 1);
end
  
Fr = Fr / max(abs(Fz(1,:)));
Ft = Ft / max(abs(Fz(1,:)));
Fz = Fz / max(abs(Fz(1,:)));

fh = fopen(strcat(fname,'.opal'),'w');
fprintf(fh, "2DDynamic XZ\n");
fprintf(fh, "% 5.10e    % 5.10e    % 6d\n", z_values(1), z_values(Nz), Nz-1);
fprintf(fh, "% 5.10e\n", FREQUENCY);
fprintf(fh, "% 5.10e    % 5.10e    % 6d\n", r_values(1), r_values(Nr), Nr-1); 
for i = 1:Nr
  for j = 1:Nz
    fprintf(fh, "% 5.10e    % 5.10e    % 5.10e    % 5.10e\n", Fz(i,j), Fr(i,j), 0.0, Ft(i,j));
  end
end

fclose(fh);
