#! /usr/bin/perl

use warnings;
use strict;

use lib "scripts/netcdf";
use abinit_readtype;

my $type = "";
my $cont = 0;
my $wait = 1;

my ($var,$grp,$i,@tmp1,@tmp2,%vars,%dims);

%vars = extract_variables("Src_defs/defs_datatypes.F90","dataset_type");

#print "Writing backup routine...";

open(OUT,">Src_3ionetcdf/def_dtset_ncdf.F90")
or die "Could not open Src_3ionetcdf/def_dtset_ncdf.F90 for writing";

print OUT <<EOF
!{\\src2tex{textfont=tt}}
!!****f* ABINIT/def_dtset_ncdf
!! NAME
!! def_dtset_ncdf
!!
!! FUNCTION
!!  Prepares a specified NetCDF file to store a dataset.
!!
!! COPYRIGHT
!!  Copyright (C) 2005 ABINIT group (YP)
!!  This file is distributed under the terms of the
!!  GNU General Public License, see ~ABINIT/Infos/copyright
!!  or http://www.gnu.org/copyleft/gpl.txt .
!!  For the initials of contributors, see ~ABINIT/Infos/contributors .
!!
!! INPUTS
!!  ncid = identifier of an open NetCDF file
!!  dims <type(ncdims_type)> = dimensions of dataset variables
!!
!! OUTPUT
!!  (no direct output : data written to a NetCDF file)
!!
!! NOTES
!!  Pouring data into the file is performed by the put_dtset_ncdf routine.
!!
!!  Please do not edit this file. Since it is automatically generated (see
!!  the makencdf script), your changes would be systematically overwritten.
!!  Instead you should directly modify the generating script itself
!!  (~ABINIT/scripts/netcdf/def_dtset_ncdf.pl) to suit your needs.
!!
!! PARENTS
!!
!! CHILDREN
!!
!! SOURCE

subroutine def_dtset_ncdf(ncid,dims)

 use defs_basis
 use defs_datatypes
#if defined NETCDF
 use netcdf
 use defs_netcdf
#endif

 implicit none

!Arguments ------------------------------------
#if defined NETCDF
 type(ncdims_type),intent(in) :: dims
#else
 integer          ,intent(in) :: dims
#endif
!scalars
 integer,intent(in) :: ncid

!Local variables-------------------------------
!scalars
 integer :: ncerr,ncvar
 character(len=500) :: message
!arrays
 integer,allocatable :: ncdims(:)

! ***************************************************************************
#if defined NETCDF
EOF
;

foreach $var (sort keys %vars)
{
	SWITCH:
	{
		( $vars{$var} eq "integer" ) && do{ $type = "NF90_INT"; last SWITCH; };
		( $vars{$var} eq "real(dp)" ) && do{ $type = "NF90_DOUBLE"; last SWITCH; };
	}

	if ( $var =~ /^([^\(]*)\(([^\)]*)\)/ )
	{
		$var  = $1;
		@tmp1 = split(/,/,$2);

		printf OUT "  ! Define dimensions of %s\n",$var;
		printf OUT "  allocate(ncdims(%d))\n",$#tmp1+1;

		for ($i=0;$i<=$#tmp1;$i++)
		{
			$_ = $tmp1[$i];
			if ( /[a-z]/ )
			{
				s/\*/_/g;
				$tmp1[$i] =~ s/\*/\*dims\%/g;
				$tmp1[$i] = "dims\%".$tmp1[$i];

				if ( $dims{$_} )
				{
					printf OUT "  ncerr = nf90_inq_dimid(ncid,name='%s',dimid=ncdims(%d))\n",$_,$i+1;
					print OUT "  call handle_ncerr_new('INQ','DIM','$_','def_dtset_ncdf',ncerr)\n";
					$dims{$_}++;
				}
				else
				{
					printf OUT "  ncerr = nf90_def_dim(ncid,name='%s',len=%s,dimid=ncdims(%d))\n",$_,$tmp1[$i],$i+1;
					print OUT "  call handle_ncerr_new('DEF','DIM','$_','def_dtset_ncdf',ncerr)\n";
					$dims{$_} = 1;
				}
			}
			else
			{
				printf OUT "  ncerr = nf90_def_dim(ncid,name='%s%d',len=%s,dimid=ncdims(%d))\n",$var."_dim",$i+1,$_,$i+1;
				printf OUT "  call handle_ncerr_new('DEF','DIM','%s%d','def_dtset_ncdf',ncerr)\n",
				 $var."_dim",$i+1;
			}
		}

		print OUT "  ! Define $var\n";
		print OUT "  ncerr = nf90_def_var(ncid,name='$var',xtype=$type,dimids=ncdims,varid=ncvar)\n";
		print OUT "  deallocate(ncdims)\n";
	}
	else
	{
		print OUT "  ! Define $var\n";
		print OUT "  ncerr = nf90_def_var(ncid,name='$var',xtype=$type,varid=ncvar)\n";
	}

	print OUT "  call handle_ncerr_new('DEF','VAR','$var','def_dtset_ncdf',ncerr)\n";
}

print OUT <<EOF
#else
 write(message,'(4a)') ch10,' def_dtset_ncdf - ERROR',ch10,&
&  'This version has not been compiled with NetCDF support'
 call wrtout(std_out,message,'COLL')
 call leave_new('COLL')
#endif

end subroutine def_dtset_ncdf
!!***
EOF
;

close(OUT);

#print "done.\n\n";

