function [gdat_data] = get_grids_1d(gdat_data,nbdim_x,nopt,nverbose);
%
% add various rhos in grids_1d
%
% assume gdat_data.x=rhopol, gdat_data.t, and gdat_data(x,t)
%
% nbdim_x = 1: same x(1:length(x)) (rhopol) for all times
%         = 2: x(:,t), rhopol depends on time like Thomson projection
%
% nopt = 0: do not fill in, just make the empty structure
%      = 1: do compute various fields
%
% compute psi, rhotor_edge, rhotornorm, volume_edge and rhovol
%
% use gdat calls to psi_axis, psi_edge, rhotor, etc with same basic parameters as data.gdat_params
%

gdat_data.grids_1d.rhopolnorm = gdat_data.x;
if (nopt == 0) || isempty(gdat_data.x) || isempty(gdat_data.t) || isempty(gdat_data.data) || ischar(gdat_data.data)
  gdat_data.grids_1d.rhotornorm = [];
  gdat_data.grids_1d.rhovolnorm = [];
  gdat_data.grids_1d.psi = [];
  gdat_data.grids_1d.rhotor_edge = [];
  gdat_data.grids_1d.volume_edge = [];
  return
end

params_eff = gdat_data.gdat_params;params_eff.doplot=0;
params_eff.data_request='rhotor_norm';
rhotor_norm = gdat(gdat_data.shot,params_eff);
ndim_x_rhotor = length(find(size(rhotor_norm.x)>1));
params_eff.data_request='rhovol';
rhovol = gdat(gdat_data.shot,params_eff);
ndim_x_rhovol = length(find(size(rhotor_norm.x)>1));
params_eff.data_request='psi_axis';
psi_axis = gdat(gdat_data.shot,params_eff);
params_eff.data_request='psi_edge';
psi_edge = gdat(gdat_data.shot,params_eff);
psi0_edge = interpos(63,psi_axis.t,psi_axis.data-psi_edge.data,gdat_data.t,-0.01);
if (nbdim_x == 1)
  gdat_data.grids_1d.psi = gdat_data.grids_1d.rhopolnorm.^2*reshape(psi0_edge,1,length(psi0_edge));
elseif (nbdim_x == 2)
  gdat_data.grids_1d.psi = gdat_data.grids_1d.rhopolnorm.^2.*repmat(reshape(psi0_edge,1,length(psi0_edge)),size(gdat_data.grids_1d.rhopolnorm,1),1);
else
  if nverbose>=0; disp(['option: nbdim_x = ' numstr(nbdim_x) ' not implemented, check with O. Sauter']); end
  return
end
gdat_data.grids_1d.rhotornorm = NaN*ones(size(gdat_data.data));
gdat_data.grids_1d.rhovolnorm = NaN*ones(size(gdat_data.data));

if (isempty(rhotor_norm.x) ||isempty(rhotor_norm.t) || isempty(rhotor_norm.data)) ...
      && (isempty(rhovol.x) ||isempty(rhovol.t) || isempty(rhovol.data))
  return
end

it_rt = iround_os(rhotor_norm.t,gdat_data.t);
it_vol = iround_os(rhovol.t,gdat_data.t);
for it=1:length(gdat_data.t)
  % do an interpolation on closest point to avoid 2D interp
  it_rt_eff = it_rt(it);
  it_vol_eff = it_vol(it);
  if (nbdim_x == 1)
    ii=find(isfinite(gdat_data.grids_1d.rhopolnorm));
  else
    ii=find(isfinite(gdat_data.grids_1d.rhopolnorm(:,it)));
  end
  if (nbdim_x == 1)
    if length(ii)==length(gdat_data.grids_1d.rhopolnorm)
      if ~isempty(rhotor_norm.x) && ~isempty(rhotor_norm.data)
        if ndim_x_rhotor==1
          gdat_data.grids_1d.rhotornorm(:,it)=interpos(-13,rhotor_norm.x,rhotor_norm.data(:,it_rt_eff),gdat_data.grids_1d.rhopolnorm);
        else
          gdat_data.grids_1d.rhotornorm(:,it)=interpos(-13,rhotor_norm.x(:,it_rt_eff),rhotor_norm.data(:,it_rt_eff),gdat_data.grids_1d.rhopolnorm);
        end
      end
      if ~isempty(rhovol.x) && ~isempty(rhovol.data)
	if ndim_x_rhovol==1
	  gdat_data.grids_1d.rhovolnorm(:,it)=interpos(-13,rhovol.x,rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm);
	else
	  gdat_data.grids_1d.rhovolnorm(:,it)=interpos(-13,rhovol.x(:,it_vol_eff),rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm);
	end
      end
    end
  else
    if length(ii)==size(gdat_data.grids_1d.rhopolnorm,1)
      if ~isempty(rhotor_norm.x) && ~isempty(rhotor_norm.data)
        if ndim_x_rhotor==1
          gdat_data.grids_1d.rhotornorm(:,it)=interpos(-13,rhotor_norm.x,rhotor_norm.data(:,it_rt_eff),gdat_data.grids_1d.rhopolnorm(:,it));
        else
          gdat_data.grids_1d.rhotornorm(:,it)=interpos(-13,rhotor_norm.x(:,it_rt_eff),rhotor_norm.data(:,it_rt_eff),gdat_data.grids_1d.rhopolnorm(:,it));
        end
      end
      if ~isempty(rhovol.x) && ~isempty(rhovol.data)
	if ndim_x_rhovol==1
	  gdat_data.grids_1d.rhovolnorm(:,it)=interpos(-13,rhovol.x,rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm(:,it));
	else
	  gdat_data.grids_1d.rhovolnorm(:,it)=interpos(-13,rhovol.x(:,it_vol_eff),rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm(:,it));
	end
      end
    end
  end
end
gdat_data.grids_1d.rhotor_edge=interpos(-63,rhotor_norm.t',rhotor_norm.rhotor_edge,gdat_data.t',-0.01);
gdat_data.grids_1d.volume_edge=interpos(-63,rhovol.t',rhovol.volume_edge,gdat_data.t',-0.01);


