problem with put_att (maximum number of attributes)

Wei-keng Liao wkliao at ece.northwestern.edu
Tue Feb 24 11:30:02 CST 2009


Lex,
I think it may be the Fortran compiler you used.
I am using gfortran, gcc version 4.1.2 and have no problem even if
I used (/istart,jstart/),(/nx_mpi,ny_mpi/) in the arguments.

Rob, can you explain why this is happening? I wonder if it is because
our Fortran interfaces could not handle such arguments.

Lex got the MPI error (see below) when he used

     nfmpi_put_vara_real(NCIDP,array_varIDP,(/istart,jstart/),(/ 
nx_mpi,ny_mpi/), buf)

but no error if he used array arguments, starts(2) and counts(2)

     integer(MPI_OFFSET_KIND) :: starts(2), counts(2)
     ....
     starts = (/istart,jstart/)
     counts = (/nx_mpi,ny_mpi/)
     nfmpi_put_vara_real_all(NCIDP,array_varIDP,starts,counts,buf)

Wei-keng


On Feb 24, 2009, at 11:15 AM, Alex Hoffmann wrote:

> Yes, this works, also if I replace in the old version the vectors for
> start and count by corresponding vector variables.  Thank you.
> Is this a function interface issue?
> Lex
>
> Wei-keng Liao wrote:
>> Hi, Lex,
>>
>> Can you try the attached program? I removed the serial netCDF part.
>> Since I don't have any problem of running it, I would suggest you to
>> run mpiexec in the debug mode and find the exact error location in  
>> pnetCDF.
>> Please let us know.
>>
>> Wei-keng
>>
>>
>>
>>
>>
>>
>> On Feb 24, 2009, at 3:33 AM, Alex Hoffmann wrote:
>>
>>> Dear Wei-keng,
>>>
>>> The following are the MPI errors I get, also if I subdivide the line
>>> (surely it can't be necessary to divide lines, though?).  I never  
>>> had
>>> any compilation problems...
>>>
>>> The job aborts at the nfmpi_put_vara_real_all() call, so my serial
>>> netCDF subfiles are still assigned and the parallel netCDF file  
>>> still
>>> created.
>>>
>>> Cheers,
>>>
>>> Lex
>>>
>>> aborting job:
>>> Fatal error in MPI_Waitall: Other MPI error, error stack:
>>> MPI_Waitall(273): MPI_Waitall(count=5, req_array=0xa11220,
>>> status_array=0xa1dc70) failed
>>> MPIDI_CH3_Progress_wait(207): an error occurred while handling an  
>>> event
>>> returned by MPIDU_Sock_Wait()
>>> MPIDI_CH3I_Progress_handle_sock_event(540): ch3|sock|immedread  
>>> 0xa07e50
>>> 0xa116e0 0xa0cd70
>>> MPIDU_Sock_readv(331): the supplied buffer contains invalid memory
>>> (set=0,sock=3,errno=14:Bad address)
>>> aborting job:
>>> Fatal error in MPI_Alltoall: Other MPI error, error stack:
>>> MPI_Alltoall(814): MPI_Alltoall(sbuf=0xa1b180, scount=1, MPI_INT,
>>> rbuf=0xa1b160, rcount=1, MPI_INT, comm=0x84000001) failed
>>> MPI_Waitall(273): MPI_Waitall(count=8, req_array=0xa1b240,
>>> status_array=0xa1b270) failed
>>> MPIDI_CH3_Progress_wait(207): an error occurred while handling an  
>>> event
>>> returned by MPIDU_Sock_Wait()
>>> MPIDI_CH3I_Progress_handle_sock_event(492):
>>> connection_recv_fail(1728):
>>> MPIDU_Socki_handle_read(590): connection closed by peer  
>>> (set=0,sock=1)
>>> aborting job:
>>> Fatal error in MPI_Allreduce: Other MPI error, error stack:
>>> MPI_Allreduce(1058): MPI_Allreduce(sbuf=0x7fff2f18140c,
>>> rbuf=0x7fff2f181408, count=1, MPI_INT, MPI_MAX, comm=0x84000000)  
>>> failed
>>> MPIR_Allreduce(317):
>>> MPIC_Sendrecv(161):
>>> MPIC_Wait(308):
>>> MPIDI_CH3_Progress_wait(207): an error occurred while handling an  
>>> event
>>> returned by MPIDU_Sock_Wait()
>>> MPIDI_CH3I_Progress_handle_sock_event(492):
>>> connection_recv_fail(1728):
>>> MPIDU_Socki_handle_read(614): connection failure
>>> (set=0,sock=2,errno=104:Connection reset by peer)
>>> rank 0 in job 37  Herzog131_50855   caused collective abort of all  
>>> ranks
>>> exit status of rank 0: return code 13
>>>
>>>
>>> Wei-keng Liao wrote:
>>>> Hi, Lex,
>>>>
>>>> I got a compile error of "Syntax error in array constructor" from  
>>>> the
>>>> line calling
>>>> nfmpi_put_vara_real_all(). But after divide this line into two  
>>>> lines:
>>>> call
>>>> handle_err_nc_p(nfmpi_put_vara_real_all(NCIDP,array_varIDP,(/ 
>>>> istart,jstart/),
>>>>
>>>> &
>>>>
>>>> (/nx_mpi,ny_mpi/),reshape(array,(/nx_mpi,ny_mpi/))),10)
>>>>
>>>> Everything compiles, runs OK.
>>>> Please tell us what MPI errors you have got.
>>>>
>>>> The followings is from the command-line output and ncmpidump:
>>>>
>>>> ----------------------------------------------------------------------------------------
>>>>
>>>>
>>>> % mpiexec -n 4 mpi_setup_netCDFout
>>>> Myrank 0:            0
>>>> MPI initialized!
>>>> Number of processes [nprocs]:            4
>>>> Coords of process            0  :           0           0
>>>> whole 4*4 array on each process:
>>>> 11.00000       12.00000       13.00000       14.00000
>>>> 21.00000       22.00000       23.00000       24.00000
>>>> 31.00000       32.00000       33.00000       34.00000
>>>> 41.00000       42.00000       43.00000       44.00000
>>>>
>>>> subgrid array on process with coords            0           0 :
>>>> shape:            2           2
>>>> 11.00000       12.00000
>>>> 21.00000       22.00000
>>>>
>>>>
>>>> array in nc_out_p subroutine
>>>> 11.00000       21.00000       12.00000       22.00000
>>>>
>>>> (/nx,ny/)                    2           2
>>>> (/istart,jstart/)                     1                    1
>>>>
>>>> pnetCDF var dimensions in xy                      
>>>> 4                    4
>>>> ntx_mpi                     4
>>>> ntx                4
>>>> nx_mpi                      2
>>>> nx                 2
>>>>
>>>>
>>>> -------------------------------------------------------------------------------------------------------
>>>>
>>>>
>>>> % ncmpidump out_glb.nc
>>>> netcdf out_glb {
>>>> dimensions:
>>>> x = 4 ;
>>>> y = 4 ;
>>>> variables:
>>>> float x(x) ;
>>>> float y(y) ;
>>>> float array(y, x) ;
>>>> data:
>>>>
>>>> x = 0, 0, 0, 0 ;
>>>>
>>>> y = 0, 0, 0, 0 ;
>>>>
>>>> array =
>>>> 11, 21, 31, 41,
>>>> 12, 22, 32, 42,
>>>> 13, 23, 33, 43,
>>>> 14, 24, 34, 44 ;
>>>> }
>>>>
>>>> --------------------------------------------------------------------------------------------------------
>>>>
>>>>
>>>>
>>>> On Feb 23, 2009, at 5:57 AM, Alex Hoffmann wrote:
>>>>
>>>>> Dear Wei-keng,
>>>>>
>>>>> below is the script section in collective data mode:
>>>>>
>>>>> -------
>>>>> ! this should have been possible in collective data mode, but did
>>>>> not seem to work !!!
>>>>> !    call handle_err_nc_p(nfmpi_begin_indep_data(NCIDP),9)
>>>>>
>>>>> !    call
>>>>> handle_err_nc_p(nfmpi_put_vara_real(NCIDP,array_varIDP,(/ 
>>>>> istart,jstart/),(/nx_mpi,ny_mpi/),reshape(array,(/ 
>>>>> nx_mpi,ny_mpi/))),10)
>>>>>
>>>>>
>>>>>
>>>>> call
>>>>> handle_err_nc_p(nfmpi_put_vara_real_all(NCIDP,array_varIDP,(/ 
>>>>> istart,jstart/),(/nx_mpi,ny_mpi/),reshape(array,(/ 
>>>>> nx_mpi,ny_mpi/))),10)
>>>>>
>>>>>
>>>>>
>>>>> !    call handle_err_nc_p(nfmpi_end_indep_data(NCIDP),11)
>>>>>
>>>>> call handle_err_nc_p(nfmpi_close(NCIDP),12)
>>>>> -------
>>>>>
>>>>> Basically, I replaced the independent call by the collective  
>>>>> call for
>>>>> all processes, and commented out the call into/out of  
>>>>> independent data
>>>>> mode.  This produces MPI error messages...
>>>>>
>>>>> Cheers,
>>>>>
>>>>> Lex
>>>>>
>>>>> Wei-keng Liao wrote:
>>>>>> Hi, Lex,
>>>>>>
>>>>>> Have you tried nfmpi_put_vara_real_all() in collective mode?
>>>>>> Collective mode is the default mode.
>>>>>>
>>>>>> If yes, can you show me how you call it?
>>>>>>
>>>>>> Wei-keng
>>>>>>
>>>>>>
>>>>>> On Feb 20, 2009, at 1:24 PM, Alex Hoffmann wrote:
>>>>>>
>>>>>>> Thank you for the clarifications.
>>>>>>>
>>>>>>> I included the mpi setup and netCDF output test code and the
>>>>>>> corresponding Makefile as attachments.  Maybe it can be of use  
>>>>>>> as a
>>>>>>> simple example.
>>>>>>>
>>>>>>> I am currently running my code with the bugs in the output file
>>>>>>> fixed,
>>>>>>> and it seems to be working just fine.
>>>>>>>
>>>>>>> Thanks again,
>>>>>>>
>>>>>>> Lex
>>>>>>>
>>>>>>> Wei-keng Liao wrote:
>>>>>>>> Hi, Lex, Please read below.
>>>>>>>>
>>>>>>>> On Feb 20, 2009, at 5:13 AM, Alex Hoffmann wrote:
>>>>>>>>
>>>>>>>>> Many thanks, Wei-keng, the issue seems to have been tracked  
>>>>>>>>> down...
>>>>>>>>>
>>>>>>>>> For reasons that are still not completely clear to me, using
>>>>>>>>> trim() for
>>>>>>>>> trimming the coordinate variable names and units (lines  
>>>>>>>>> 480/81,
>>>>>>>>> 485/86)
>>>>>>>>> seemed to leave a strange character (such as y diacritic or Y
>>>>>>>>> accent)
>>>>>>>>> trailing on e.g. the units (such as degrees_east).  Though  
>>>>>>>>> this
>>>>>>>>> did not
>>>>>>>>> affect the run on its own, any subsequent attribute definition
>>>>>>>>> made the
>>>>>>>>> code crash.  I fixed this by removing the trim().
>>>>>>>>
>>>>>>>> If you called trim() on a string argument in a function, the  
>>>>>>>> tailing
>>>>>>>> character '\0' will not be passed to the function. So if the  
>>>>>>>> string
>>>>>>>> length
>>>>>>>> argument has a larger value, the function will include  
>>>>>>>> undefined
>>>>>>>> tailing
>>>>>>>> characters as part of the string. Glad to know you solve the
>>>>>>>> problem.
>>>>>>>>
>>>>>>>>
>>>>>>>>> If I may, I'd like to clarify one further issue with PnetCDF.
>>>>>>>>> I initially implemented the code (also a small test program)  
>>>>>>>>> in
>>>>>>>>> collective data mode.  This did not work (the different
>>>>>>>>> processes/Cartesian sub domains) did not write their data  
>>>>>>>>> into the
>>>>>>>>> netCDF file.  I then changed the data write statements into
>>>>>>>>> independent
>>>>>>>>> data mode, with which it worked immediately.  Am I right to  
>>>>>>>>> use
>>>>>>>>> independent data mode for this purpose?
>>>>>>>>
>>>>>>>> I am not sure what you meant here.
>>>>>>>> Can you post your small test program to the list?
>>>>>>>>
>>>>>>>> Wei-keng
>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> Thank you very much for your support.
>>>>>>>>>
>>>>>>>>> Cheers
>>>>>>>>>
>>>>>>>>> Lex
>>>>>>>>>
>>>>>>>>> Wei-keng Liao wrote:
>>>>>>>>>> Hi, Lex,
>>>>>>>>>>
>>>>>>>>>> If this error message "NC definations on multiprocesses  
>>>>>>>>>> conflict"
>>>>>>>>>> is the only one you got before every crash, it is most  
>>>>>>>>>> likely the
>>>>>>>>>> attribute values are not the same across all processes.
>>>>>>>>>> Please check that first, so you can rule out that  
>>>>>>>>>> possibility.
>>>>>>>>>>
>>>>>>>>>> Also, can you check the attributes in the output file for  
>>>>>>>>>> the case
>>>>>>>>>> that ran successfully? (using ncmpidump -c filename) This  
>>>>>>>>>> is just
>>>>>>>>>> to see if the attributes are the ones you expected. Also,  
>>>>>>>>>> please
>>>>>>>>>> make sure the number of attributes created is less than 4K  
>>>>>>>>>> global
>>>>>>>>>> or per variable, as Rob indicated.
>>>>>>>>>>
>>>>>>>>>> Wei-keng
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> On Feb 19, 2009, at 5:53 PM, Alex Hoffmann wrote:
>>>>>>>>>>
>>>>>>>>>>> Dear Wei-keng,
>>>>>>>>>>>
>>>>>>>>>>> I understand that the attribute define arguments must be  
>>>>>>>>>>> the same
>>>>>>>>>>> across
>>>>>>>>>>> all processes.   I am also almost certain that the  
>>>>>>>>>>> arguments you
>>>>>>>>>>> cited
>>>>>>>>>>> are identical across processes but I will perform a check.
>>>>>>>>>>>
>>>>>>>>>>> Be this as it may, I get the same arrow message if I  
>>>>>>>>>>> uncomment
>>>>>>>>>>> one of
>>>>>>>>>>> the text attributes specifying the unit or one of the real
>>>>>>>>>>> attributes
>>>>>>>>>>> specifying the fill value, for instance.  There is no  
>>>>>>>>>>> possibility
>>>>>>>>>>> that
>>>>>>>>>>> these are different across processes as they are introduced
>>>>>>>>>>> straight
>>>>>>>>>>> into the function call.  Any check I performed seemed to  
>>>>>>>>>>> indicate
>>>>>>>>>>> that
>>>>>>>>>>> the problem arises from the number of arguments I defined,  
>>>>>>>>>>> and
>>>>>>>>>>> not
>>>>>>>>>>> from
>>>>>>>>>>> the arguments themselves, since they work all fine  
>>>>>>>>>>> individually.
>>>>>>>>>>>
>>>>>>>>>>> I was actually thinking to write the attributes ad-hoc  
>>>>>>>>>>> once the
>>>>>>>>>>> dataset
>>>>>>>>>>> is defined using regular netCDF functions.  This however
>>>>>>>>>>> seemed to
>>>>>>>>>>> be a
>>>>>>>>>>> little elegant solution defeating the purpose of using  
>>>>>>>>>>> PnetCDF.
>>>>>>>>>>>
>>>>>>>>>>> Thanks,
>>>>>>>>>>>
>>>>>>>>>>> Lex
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> Wei-keng Liao a écrit :
>>>>>>>>>>>> Hi, Lex,
>>>>>>>>>>>>
>>>>>>>>>>>> PnetCDF requires the attribute define arguments be the same
>>>>>>>>>>>> across
>>>>>>>>>>>> all processes. This is why you got the error message for
>>>>>>>>>>>> defining
>>>>>>>>>>>> attribute "history" in line 400, if its value "conflicts"  
>>>>>>>>>>>> among
>>>>>>>>>>>> processes.
>>>>>>>>>>>> I think you have got this.
>>>>>>>>>>>>
>>>>>>>>>>>> I checked your codes. Starting from line 562, you have a  
>>>>>>>>>>>> few
>>>>>>>>>>>> text
>>>>>>>>>>>> attributes
>>>>>>>>>>>> from des_tgas_picture(i), des_trac_picture(i),
>>>>>>>>>>>> des_tpas_picture(i),
>>>>>>>>>>>> and des.
>>>>>>>>>>>> Similar for out_type == 2. Are you sure their values are  
>>>>>>>>>>>> all the
>>>>>>>>>>>> same
>>>>>>>>>>>> across
>>>>>>>>>>>> all the processes? May be you can use MPI_Bcast to check if
>>>>>>>>>>>> anyone has
>>>>>>>>>>>> different value.
>>>>>>>>>>>>
>>>>>>>>>>>> Wei-keng
>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>> On Feb 19, 2009, at 11:14 AM, Alex Hoffmann wrote:
>>>>>>>>>>>>
>>>>>>>>>>>>> Dear Rob
>>>>>>>>>>>>> Thanks a lot for your answer.
>>>>>>>>>>>>> Isolating the PnetCDF output module from the code as a
>>>>>>>>>>>>> stand-alone
>>>>>>>>>>>>> would
>>>>>>>>>>>>> require a substantial effort, and I do not encounter  
>>>>>>>>>>>>> this error
>>>>>>>>>>>>> in  my
>>>>>>>>>>>>> small test programs.  I attached the output module for
>>>>>>>>>>>>> reference.
>>>>>>>>>>>>>
>>>>>>>>>>>>> The attribute section that causes trouble is between lines
>>>>>>>>>>>>> 369-405
>>>>>>>>>>>>> (global) and 468-613 (variables).  I am aware that the
>>>>>>>>>>>>> attribute
>>>>>>>>>>>>> history
>>>>>>>>>>>>> (line 400) may be problematic as every process may have a
>>>>>>>>>>>>> slightly
>>>>>>>>>>>>> different value (for the time element taken from the  
>>>>>>>>>>>>> fortran
>>>>>>>>>>>>> function
>>>>>>>>>>>>> date_and_time).  Yet this is not the only issue.
>>>>>>>>>>>>>
>>>>>>>>>>>>> If I prematurely exit at line 492, uncommenting any  
>>>>>>>>>>>>> further
>>>>>>>>>>>>> attribute
>>>>>>>>>>>>> definition between 468 and 489 causes a crash with the
>>>>>>>>>>>>> multiprocesses
>>>>>>>>>>>>> conflict.
>>>>>>>>>>>>>
>>>>>>>>>>>>> Note also that the attribute definitions between 558 and  
>>>>>>>>>>>>> 613
>>>>>>>>>>>>> only
>>>>>>>>>>>>> produce no crash currently because of the print  
>>>>>>>>>>>>> statements on
>>>>>>>>>>>>> lines 560
>>>>>>>>>>>>> and 591, which made me think that maybe there is an
>>>>>>>>>>>>> optimization
>>>>>>>>>>>>> problem.  What is strange is that if these are not  
>>>>>>>>>>>>> included,
>>>>>>>>>>>>> the
>>>>>>>>>>>>> crash
>>>>>>>>>>>>> occurs only when I start writing the first file (movie  
>>>>>>>>>>>>> file,
>>>>>>>>>>>>> out_type=2), AND if I leave the first block (if  
>>>>>>>>>>>>> out_type==1)
>>>>>>>>>>>>> uncommented, although this block is not even entered in  
>>>>>>>>>>>>> that
>>>>>>>>>>>>> particular
>>>>>>>>>>>>> case, and vice versa.  Hope this is more or less clear...
>>>>>>>>>>>>>
>>>>>>>>>>>>> Cheers,
>>>>>>>>>>>>> Lex
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> Robert Latham wrote:
>>>>>>>>>>>>>> On Thu, Feb 19, 2009 at 10:09:37AM +0000, Alex Hoffmann  
>>>>>>>>>>>>>> wrote:
>>>>>>>>>>>>>>> Hello,
>>>>>>>>>>>>>>> I am getting the same error message as had already been
>>>>>>>>>>>>>>> inquired
>>>>>>>>>>>>>>> on in
>>>>>>>>>>>>>>> this list, i.e.:
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> NC definations on multiprocesses conflict.
>>>>>>>>>>>>>>> NOTE: Definitions across all process
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> and this seems to be related somehow randomly to the  
>>>>>>>>>>>>>>> total
>>>>>>>>>>>>>>> number of
>>>>>>>>>>>>>>> attributes (global or for each variable) that I intend  
>>>>>>>>>>>>>>> to
>>>>>>>>>>>>>>> write
>>>>>>>>>>>>>>> during
>>>>>>>>>>>>>>> the the definition phase (and, though I am not sure  
>>>>>>>>>>>>>>> about
>>>>>>>>>>>>>>> this, to
>>>>>>>>>>>>>>> the
>>>>>>>>>>>>>>> number of processes used to run my F90 parallel code).
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> More precisely, if I uncomment one put_att to replace
>>>>>>>>>>>>>>> another,
>>>>>>>>>>>>>>> the
>>>>>>>>>>>>>>> program still runs, if I do not comment out another  
>>>>>>>>>>>>>>> one, the
>>>>>>>>>>>>>>> program
>>>>>>>>>>>>>>> crashes with the above error message.
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Is there any way of identifying where this problem  
>>>>>>>>>>>>>>> comes from
>>>>>>>>>>>>>>> or is
>>>>>>>>>>>>>>> there a maximum amount of attributes (though I don't  
>>>>>>>>>>>>>>> write an
>>>>>>>>>>>>>>> excessive
>>>>>>>>>>>>>>> amount of them)?  There is no problem running it on one
>>>>>>>>>>>>>>> process(or).
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Maybe you can show us with a small test program how you  
>>>>>>>>>>>>>> are
>>>>>>>>>>>>>> encountering this error?  The maximum number of  
>>>>>>>>>>>>>> attributes
>>>>>>>>>>>>>> is 4k
>>>>>>>>>>>>>> (the
>>>>>>>>>>>>>> constant NC_MAX_ATTRS).
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ==rob
>>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> -- 
>>>>>>>>>>>>> ______________________________________________________________________
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> Alex HOFFMANN                     PhD Candidate
>>>>>>>>>>>>> Centre for Atmospheric Science
>>>>>>>>>>>>> Department of Geography           University of Cambridge
>>>>>>>>>>>>> Downing Place, Cambridge CB2 3EN, UK
>>>>>>>>>>>>> e-mail: ah519 at cam.ac.uk  tel: +44 (0)1223 766581
>>>>>>>>>>>>> www.geog.cam.ac.uk/people/hoffmann/
>>>>>>>>>>>>> ______________________________________________________________________
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> !-*- F90 -*-
>>>>>>>>>>>>> module atham_pnetcdf
>>>>>>>>>>>>>
>>>>>>>>>>>>> !--------------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! author:  Alex
>>>>>>>>>>>>> HOFFMANN                                             !
>>>>>>>>>>>>> ! email:
>>>>>>>>>>>>> ah519 
>>>>>>>>>>>>> @cam.ac.uk                                           !
>>>>>>>>>>>>> ! date:    February
>>>>>>>>>>>>> 2009                                             !
>>>>>>>>>>>>> ! version:
>>>>>>>>>>>>> v0.1 
>>>>>>>>>>>>>                                                       !
>>>>>>>>>>>>> !
>>>>>>>>>>>>> !
>>>>>>>>>>>>> ! collection of routines for ATHAM PnetCDF
>>>>>>>>>>>>> output                    !
>>>>>>>>>>>>> !--------------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! module is to replace existing atham_netcdf.mod (GL) to  
>>>>>>>>>>>>> write
>>>>>>>>>>>>> ATHAM
>>>>>>>>>>>>> output as netCDF files
>>>>>>>>>>>>> ! same as atham_netcdf.mod (AH) but adapted to run with  
>>>>>>>>>>>>> the
>>>>>>>>>>>>> Parallel
>>>>>>>>>>>>> NetCDF API,
>>>>>>>>>>>>> !  developed by Northwestern University and Argonne  
>>>>>>>>>>>>> National
>>>>>>>>>>>>> Laboratory,
>>>>>>>>>>>>> !  which is built on top of MPI-IO (tested with MPICH 2.1)
>>>>>>>>>>>>>
>>>>>>>>>>>>> use precision,    only : kint, kreal, real4, my_real
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only : nx, ny, nz, ntx, nty, npx, npy
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only : ntgas_picture, ntrac_picture,
>>>>>>>>>>>>> ntpas_picture, &
>>>>>>>>>>>>>                 ntgas_movie, ntrac_movie, ntpas_movie
>>>>>>>>>>>>>
>>>>>>>>>>>>> implicit none
>>>>>>>>>>>>>
>>>>>>>>>>>>> private
>>>>>>>>>>>>> public :: netcdf_data
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! 'pnetcdf.inc' file includes cpp directives that were not
>>>>>>>>>>>>> understood
>>>>>>>>>>>>> at compile-time
>>>>>>>>>>>>> ! this is to be changed and hard-wired link to be replaced
>>>>>>>>>>>>> by the
>>>>>>>>>>>>> include 'pnetcdf.inc' command
>>>>>>>>>>>>> !include 'pnetcdf.inc'
>>>>>>>>>>>>> #include "/usr/local/parallel-netcdf-1.0.1/include/ 
>>>>>>>>>>>>> pnetcdf.inc"
>>>>>>>>>>>>> include 'mpif.h'
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! ######### tried saving netCDF variable IDs here so  
>>>>>>>>>>>>> that they
>>>>>>>>>>>>> must
>>>>>>>>>>>>> not be retrieved by inquiry
>>>>>>>>>>>>> ! ######### if this creates conflict, especially between
>>>>>>>>>>>>> picture
>>>>>>>>>>>>> and
>>>>>>>>>>>>> movie file, this must be changed
>>>>>>>>>>>>> integer(kint), save                   :: x_dimID, y_dimID,
>>>>>>>>>>>>> z_dimID, t_dimID
>>>>>>>>>>>>> integer(kint), save                   :: x_varID, y_varID,
>>>>>>>>>>>>> z_varID, t_varID
>>>>>>>>>>>>> integer(kint), save                   :: u_varID, v_varID,
>>>>>>>>>>>>> w_varID
>>>>>>>>>>>>> integer(kint), save                   :: pnew_varID,
>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>> tempnew_varID, density_varID
>>>>>>>>>>>>> integer(kint), save                   :: turbhor_varID,
>>>>>>>>>>>>> turbver_varID, turblen_varID
>>>>>>>>>>>>> integer(kint), dimension(99), save    :: tgas_varID,
>>>>>>>>>>>>> trac_varID,
>>>>>>>>>>>>> tpas_varID, d_trac_varID
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> integer(MPI_OFFSET_KIND), save        :: tcounter
>>>>>>>>>>>>>
>>>>>>>>>>>>> real(real4), parameter                :: undef =
>>>>>>>>>>>>> -0.99e33_real4 !
>>>>>>>>>>>>> make sure precision of _FillValue and output  
>>>>>>>>>>>>> (NF90_FLOAT) is
>>>>>>>>>>>>> the
>>>>>>>>>>>>> same
>>>>>>>>>>>>> !!!
>>>>>>>>>>>>>
>>>>>>>>>>>>> contains
>>>>>>>>>>>>> ! 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> ==========================================================
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> subroutine netcdf_data(outfile,out_type,timetot)
>>>>>>>>>>>>> ! IF: outfile = picture_file or movie_file (path/ 
>>>>>>>>>>>>> filename for
>>>>>>>>>>>>> output created in atham_module)
>>>>>>>>>>>>> !     out_type = 1 (picture_file) or 2 (movie_file)
>>>>>>>>>>>>> !     timetot = total simulation time since simulation  
>>>>>>>>>>>>> start
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only : unew, vnew, wnew, pnew, tetnew,
>>>>>>>>>>>>> tempnew, &
>>>>>>>>>>>>>                     density, turbhor, turbver, turblen
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only : tgasnew, tracnew, tpasnew,  
>>>>>>>>>>>>> tracdep
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only : itgas_picture, itrac_picture,
>>>>>>>>>>>>> itpas_picture, &
>>>>>>>>>>>>>                     itgas_movie, itrac_movie, itpas_movie
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only : mycoor, my_cart, myrank
>>>>>>>>>>>>>
>>>>>>>>>>>>> character(*),           intent(in   ) :: outfile
>>>>>>>>>>>>> integer(kint),          intent(in   ) :: out_type
>>>>>>>>>>>>> real(kreal),            intent(in   ) :: timetot
>>>>>>>>>>>>>
>>>>>>>>>>>>> !------------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! local
>>>>>>>>>>>>> variables 
>>>>>>>>>>>>>                                                   !
>>>>>>>>>>>>> !------------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> integer(MPI_OFFSET_KIND)              :: nx_mpi,ny_mpi,  
>>>>>>>>>>>>> nz_mpi,
>>>>>>>>>>>>> y_cl_mpi
>>>>>>>>>>>>> integer(MPI_OFFSET_KIND)              :: istart,jstart
>>>>>>>>>>>>> integer(kint)                         :: STATUS
>>>>>>>>>>>>> character(len=50)                     :: OP  ! debugging
>>>>>>>>>>>>> variable
>>>>>>>>>>>>> corresponding to the called netCDF OPeration
>>>>>>>>>>>>> integer(kint)                         :: NCID
>>>>>>>>>>>>> integer(kint)                         :: i
>>>>>>>>>>>>> logical                               :: DynMovie_flag,
>>>>>>>>>>>>> bounds_control
>>>>>>>>>>>>> integer(kint)                         :: y_st, y_sp, y_cl
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! replace DynMovie precompiler flag by program flag for  
>>>>>>>>>>>>> more
>>>>>>>>>>>>> convenient code structure
>>>>>>>>>>>>> DynMovie_flag = .false.
>>>>>>>>>>>>> #ifdef DynMovie
>>>>>>>>>>>>> DynMovie_flag = .true.
>>>>>>>>>>>>> #endif
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! re-assign dimension lengths as an MPI_Offset type  
>>>>>>>>>>>>> (necessary
>>>>>>>>>>>>> for
>>>>>>>>>>>>> the PnetCDF interface)
>>>>>>>>>>>>> ! cut-off subdomain borders
>>>>>>>>>>>>> nx_mpi  = nx-2
>>>>>>>>>>>>> ny_mpi  = ny-2
>>>>>>>>>>>>> nz_mpi  = nz
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! assign coordinates of first element of each process in  
>>>>>>>>>>>>> the
>>>>>>>>>>>>> netCDF file
>>>>>>>>>>>>> istart = mycoor(1)*nx_mpi + 1_MPI_OFFSET_KIND
>>>>>>>>>>>>> jstart = mycoor(2)*ny_mpi + 1_MPI_OFFSET_KIND
>>>>>>>>>>>>>
>>>>>>>>>>>>> bounds_control = .true. ! debugging flag
>>>>>>>>>>>>> tcounter = tcounter+1_MPI_OFFSET_KIND
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! === open netCDF dataset: enter data mode
>>>>>>>>>>>>> ------------------------!
>>>>>>>>>>>>> ! if     dataset does not yet exist (STATUS =  ), create
>>>>>>>>>>>>> file and
>>>>>>>>>>>>> enter data mode
>>>>>>>>>>>>> ! else   continue
>>>>>>>>>>>>> OP = 'open existing netCDF file'
>>>>>>>>>>>>> STATUS = nfmpi_open(my_cart, trim(outfile), NF_WRITE,
>>>>>>>>>>>>> MPI_INFO_NULL, NCID)
>>>>>>>>>>>>> if ( STATUS /= nf_noerr ) then
>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>>  print *, ""
>>>>>>>>>>>>>  print *, "P-NetCDF operation: ",trim(OP)
>>>>>>>>>>>>>  print *, trim(nfmpi_strerror(STATUS))
>>>>>>>>>>>>>  print *, "atham_pnetcdf netcdf_data(): new netCDF file
>>>>>>>>>>>>> created: ", outfile
>>>>>>>>>>>>>  print *, ""
>>>>>>>>>>>>> end if
>>>>>>>>>>>>> call netcdf_define(outfile,out_type,timetot,NCID)
>>>>>>>>>>>>> tcounter = 1_MPI_OFFSET_KIND ! reset time counter for  
>>>>>>>>>>>>> output if
>>>>>>>>>>>>> new file is created
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! ---                        write data values
>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>> ! unlimited record dimension time (only one process  
>>>>>>>>>>>>> writes the
>>>>>>>>>>>>> whole dimension, i.e. accross all processes)
>>>>>>>>>>>>> OP = 'write data values - unlimited dimension t'
>>>>>>>>>>>>>
>>>>>>>>>>>>> call handle_err(nfmpi_begin_indep_data(NCID))
>>>>>>>>>>>>>
>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID, t_varID,
>>>>>>>>>>>>> (/tcounter/), (/1_MPI_OFFSET_KIND/),  
>>>>>>>>>>>>> real(timetot,real4)),OP)
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> call handle_err(nfmpi_end_indep_data(NCID))
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! ---                        write data values
>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>> ! discriminate between 2D and 3D case for output  
>>>>>>>>>>>>> generation
>>>>>>>>>>>>> if (nty>4) then
>>>>>>>>>>>>> y_st = 2    ! y index start
>>>>>>>>>>>>> y_sp = ny-1 ! y index stop
>>>>>>>>>>>>> y_cl = ny-2 ! y number of columns
>>>>>>>>>>>>> else
>>>>>>>>>>>>> y_st = ny-1 ! 3 in 2D case
>>>>>>>>>>>>> y_sp = ny-1 ! 3 in 2D case
>>>>>>>>>>>>> y_cl = 1
>>>>>>>>>>>>> end if
>>>>>>>>>>>>> y_cl_mpi = y_cl
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! dynamics data
>>>>>>>>>>>>> OP = 'write data values - data'
>>>>>>>>>>>>> if (out_type == 1 .OR. DynMovie_flag == .true.) then
>>>>>>>>>>>>> call write_data(real(unew,real4),u_varID)
>>>>>>>>>>>>> bounds_control = .false.
>>>>>>>>>>>>> call write_data(real(vnew,   real4),v_varID)
>>>>>>>>>>>>> call write_data(real(wnew,   real4),w_varID)
>>>>>>>>>>>>> call write_data(real(pnew,   real4),pnew_varID)
>>>>>>>>>>>>> call write_data(real(tetnew, real4),tetnew_varID)
>>>>>>>>>>>>> call write_data(real(tempnew,real4),tempnew_varID)
>>>>>>>>>>>>> call write_data(real(density,real4),density_varID)
>>>>>>>>>>>>> call write_data(real(turbhor,real4),turbhor_varID)
>>>>>>>>>>>>> call write_data(real(turbver,real4),turbver_varID)
>>>>>>>>>>>>> call write_data(real(turblen,real4),turblen_varID)
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! tracer data
>>>>>>>>>>>>> if (out_type == 1) then
>>>>>>>>>>>>> do i=1,ntgas_picture
>>>>>>>>>>>>>  call
>>>>>>>>>>>>> write_data 
>>>>>>>>>>>>> (real 
>>>>>>>>>>>>> (tgasnew(:,:,:,itgas_picture(i)),real4),tgas_varID(i))
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>  call
>>>>>>>>>>>>> write_data 
>>>>>>>>>>>>> (real 
>>>>>>>>>>>>> (tracnew(:,:,:,itrac_picture(i)),real4),trac_varID(i))
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntpas_picture
>>>>>>>>>>>>>  call
>>>>>>>>>>>>> write_data 
>>>>>>>>>>>>> (real 
>>>>>>>>>>>>> (tpasnew(:,:,:,itpas_picture(i)),real4),tpas_varID(i))
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>  call handle_err(nfmpi_begin_indep_data(NCID))
>>>>>>>>>>>>>  call handle_err(nfmpi_put_vara_real(NCID,
>>>>>>>>>>>>> d_trac_varID(i),(/istart,jstart,tcounter/), &
>>>>>>>>>>>>>       (/nx_mpi,y_cl_mpi,1_MPI_OFFSET_KIND/),
>>>>>>>>>>>>> real 
>>>>>>>>>>>>> (reshape 
>>>>>>>>>>>>> (tracdep(2:nx-1,y_st:y_sp,itrac_picture(i)),shape=(/ 
>>>>>>>>>>>>> nx-2,y_cl,1/)),real4)
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> ),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_end_indep_data(NCID))
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> else if (out_type == 2) then
>>>>>>>>>>>>> do i=1,ntgas_movie
>>>>>>>>>>>>>  call
>>>>>>>>>>>>> write_data 
>>>>>>>>>>>>> (real(tgasnew(:,:,:,itgas_movie(i)),real4),tgas_varID(i))
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntrac_movie
>>>>>>>>>>>>>  call
>>>>>>>>>>>>> write_data 
>>>>>>>>>>>>> (real(tracnew(:,:,:,itrac_movie(i)),real4),trac_varID(i))
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntpas_movie
>>>>>>>>>>>>>  call
>>>>>>>>>>>>> write_data 
>>>>>>>>>>>>> (real(tpasnew(:,:,:,itpas_movie(i)),real4),tpas_varID(i))
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> endif
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! === close netCDF dataset: exit data mode & save
>>>>>>>>>>>>> -----------------!
>>>>>>>>>>>>> OP = 'close netCDF file'
>>>>>>>>>>>>> call handle_err(nfmpi_close(NCID))
>>>>>>>>>>>>>
>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>> print *, ''
>>>>>>>>>>>>> print *, 'wrote netCDF output to file ', outfile
>>>>>>>>>>>>> print *, ''
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> contains
>>>>>>>>>>>>> !----------------------------------------------------------------
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> subroutine write_data(array_in,varID)
>>>>>>>>>>>>> ! IF: array_in = data array (3D) to be written to netCDF
>>>>>>>>>>>>> !     varID    = data variable netCDF ID
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only : iflgs ! topography mask:  
>>>>>>>>>>>>> topography
>>>>>>>>>>>>> (iflgs=0) remains undef
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only  : myrank
>>>>>>>>>>>>>
>>>>>>>>>>>>> real(real4),dimension(nx,ny,nz), intent (in   )   ::  
>>>>>>>>>>>>> array_in
>>>>>>>>>>>>> integer(kint),                intent(in   )       :: varID
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! local
>>>>>>>>>>>>> variables                                              !
>>>>>>>>>>>>>
>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> real(real4), dimension(nx,ny,nz)            ::  
>>>>>>>>>>>>> array_masked
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! dimensions & bounds control
>>>>>>>>>>>>> integer(MPI_OFFSET_KIND)                    :: dimlenx,
>>>>>>>>>>>>> dimleny, dimlenz, dimlent
>>>>>>>>>>>>> integer(kint), dimension(NF_max_var_dims)   :: dimIDs
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! dimensions & bounds control
>>>>>>>>>>>>> if (myrank == 0) then
>>>>>>>>>>>>>   if (bounds_control) then
>>>>>>>>>>>>>      print *, ''
>>>>>>>>>>>>>      print *, 'atham var in xyz
>>>>>>>>>>>>> ',shape(array_in)
>>>>>>>>>>>>>      print *, 'atham var reshaped to xyzt
>>>>>>>>>>>>> ',shape(reshape(array_in(2:nx-1,y_st:y_sp,1:nz),shape=(/ 
>>>>>>>>>>>>> nx-2,y_cl,nz,1/)))
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>      call handle_err(nfmpi_inq_vardimid(NCID, u_varID,
>>>>>>>>>>>>> dimIDs)) ! taken unew as ref
>>>>>>>>>>>>>      call handle_err(nfmpi_inq_dimlen(NCID, dimIDs(1),
>>>>>>>>>>>>> dimlenx))
>>>>>>>>>>>>>      call handle_err(nfmpi_inq_dimlen(NCID, dimIDs(2),
>>>>>>>>>>>>> dimleny))
>>>>>>>>>>>>>      call handle_err(nfmpi_inq_dimlen(NCID, dimIDs(3),
>>>>>>>>>>>>> dimlenz))
>>>>>>>>>>>>>      call handle_err(nfmpi_inq_dimlen(NCID, dimIDs(4),
>>>>>>>>>>>>> dimlent))
>>>>>>>>>>>>>      print *, 'netCDF var dimensions in xyzt
>>>>>>>>>>>>> ',dimlenx,dimleny,dimlenz,dimlent
>>>>>>>>>>>>>   end if
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! mask topography
>>>>>>>>>>>>> where (iflgs == 1)
>>>>>>>>>>>>>   array_masked = array_in
>>>>>>>>>>>>> elsewhere
>>>>>>>>>>>>>   array_masked = undef
>>>>>>>>>>>>> endwhere
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! mask control
>>>>>>>>>>>>> !print *, ''
>>>>>>>>>>>>> !print *, 'topography mask'
>>>>>>>>>>>>> !print *, iflgs(1:5,3,3)
>>>>>>>>>>>>> !print *, iflgs(1:5,3,2)
>>>>>>>>>>>>> !print *, iflgs(1:5,3,1)
>>>>>>>>>>>>> !print *, ''
>>>>>>>>>>>>> !print *, 'masked variable array'
>>>>>>>>>>>>> !print *, array_masked(:,3,3)
>>>>>>>>>>>>> !print *, array_masked(:,3,2)
>>>>>>>>>>>>> !print *, array_masked(:,3,1)
>>>>>>>>>>>>> !print *, ''
>>>>>>>>>>>>>
>>>>>>>>>>>>> call handle_err(nfmpi_begin_indep_data(NCID))
>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID,
>>>>>>>>>>>>> varID,(/istart,jstart,1_MPI_OFFSET_KIND,tcounter/), &
>>>>>>>>>>>>>     (/nx_mpi,y_cl_mpi,nz_mpi,1_MPI_OFFSET_KIND/),
>>>>>>>>>>>>> reshape(array_masked(2:nx-1,y_st:y_sp,1:nz),shape=(/ 
>>>>>>>>>>>>> nx-2,y_cl,nz,1/))
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> ),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_end_indep_data(NCID))
>>>>>>>>>>>>>
>>>>>>>>>>>>> end subroutine write_data
>>>>>>>>>>>>> !
>>>>>>>>>>>>> ----------------------------------------------------------------
>>>>>>>>>>>>>
>>>>>>>>>>>>> end subroutine netcdf_data
>>>>>>>>>>>>> ! 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> ==========================================================
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> subroutine  
>>>>>>>>>>>>> netcdf_define(outfile,out_type,timetot,ncid_def)
>>>>>>>>>>>>> ! IF: outfile = picture_file or movie_file (path/ 
>>>>>>>>>>>>> filename for
>>>>>>>>>>>>> output created in atham_module)
>>>>>>>>>>>>> !     out_type = 1 (picture_file) or 2 (movie_file)
>>>>>>>>>>>>> !     timetot = simulation time since simulation start  
>>>>>>>>>>>>> (from
>>>>>>>>>>>>> atham/atham_out)
>>>>>>>>>>>>> !     ncid_def = netCDF ID assigned to file during  
>>>>>>>>>>>>> creation
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only   : var_tgas_picture,
>>>>>>>>>>>>> var_tgas_movie,           &
>>>>>>>>>>>>>                       des_tgas_picture,
>>>>>>>>>>>>> des_tgas_movie,           &
>>>>>>>>>>>>>                       var_trac_picture,
>>>>>>>>>>>>> var_trac_movie,           &
>>>>>>>>>>>>>                       des_trac_picture,
>>>>>>>>>>>>> des_trac_movie,           &
>>>>>>>>>>>>>                       var_tpas_picture,
>>>>>>>>>>>>> var_tpas_movie,           &
>>>>>>>>>>>>>                       des_tpas_picture, des_tpas_movie
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only   : iindex, spinup, nrep, periodt
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only   : volcano_setup,
>>>>>>>>>>>>> coignimbrite_setup,          &
>>>>>>>>>>>>>                       nuclear_setup, convection_setup,
>>>>>>>>>>>>> procsconfig_setup
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only   : x, y, z, xg, yg
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only   : grid2geo, sim_d_YYYY, sim_d_MM,
>>>>>>>>>>>>> sim_d_DD,   &
>>>>>>>>>>>>>                       sim_t_hh, sim_t_mm, sim_t_ss
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only   : my_cart, myrank
>>>>>>>>>>>>>
>>>>>>>>>>>>> character(*),           intent(in   ) :: outfile
>>>>>>>>>>>>> integer(kint),          intent(in   ) :: out_type
>>>>>>>>>>>>> real(kreal),            intent(in   ) :: timetot
>>>>>>>>>>>>> integer(kint),          intent(  out) :: ncid_def
>>>>>>>>>>>>>
>>>>>>>>>>>>> !------------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! local
>>>>>>>>>>>>> variables 
>>>>>>>>>>>>>                                                   !
>>>>>>>>>>>>> !------------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> integer(MPI_OFFSET_KIND)              ::  
>>>>>>>>>>>>> ntx_mpi,nty_mpi,nz_mpi
>>>>>>>>>>>>> character(len=50)                     :: OP  ! debugging
>>>>>>>>>>>>> variable
>>>>>>>>>>>>> corresponding to the called netCDF OPeration
>>>>>>>>>>>>> integer(kint)                         :: NCID
>>>>>>>>>>>>> integer(kint)                         :: i
>>>>>>>>>>>>> character(len=50)                     :: var, des
>>>>>>>>>>>>> real(kreal)                           :: toffset
>>>>>>>>>>>>> real(real4)                           :: coord_scale,
>>>>>>>>>>>>> coord_x_offset, coord_y_offset, coord_t_offset
>>>>>>>>>>>>> logical                               :: DynMovie_flag
>>>>>>>>>>>>> character(len=3)                      :: coord_x_name,
>>>>>>>>>>>>> coord_y_name, coord_z_name
>>>>>>>>>>>>> character(len=4)                      :: coord_t_name
>>>>>>>>>>>>> character(len=9)                      :: coord_x_lname,
>>>>>>>>>>>>> coord_y_lname
>>>>>>>>>>>>> character(len=13)                     :: coord_x_units,
>>>>>>>>>>>>> coord_y_units
>>>>>>>>>>>>> character(len=33)                     :: coord_t_units
>>>>>>>>>>>>> character(len=8)                      :: date
>>>>>>>>>>>>> character(len=10)                     :: time
>>>>>>>>>>>>> character(len=19)                     :: simdatetime
>>>>>>>>>>>>> character(len=10)                     :: simdate
>>>>>>>>>>>>> character(len=20)                     :: login
>>>>>>>>>>>>> character(len=63)                     :: hist_string
>>>>>>>>>>>>> integer(kint)                         :: isetup
>>>>>>>>>>>>> character(len=12)                     :: model_setup(6)
>>>>>>>>>>>>> data model_setup /'volcano     ', 'coignimbrite',
>>>>>>>>>>>>> 'nuclear     ',
>>>>>>>>>>>>> 'convection  ', 'procsconfig ','unknown     '/
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! replace DynMovie precompiler flag by program flag for  
>>>>>>>>>>>>> more
>>>>>>>>>>>>> convenient code structure
>>>>>>>>>>>>> DynMovie_flag = .false.
>>>>>>>>>>>>> #ifdef DynMovie
>>>>>>>>>>>>> DynMovie_flag = .true.
>>>>>>>>>>>>> #endif
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! re-assign dimension lengths as an MPI_Offset type  
>>>>>>>>>>>>> (necessary
>>>>>>>>>>>>> for
>>>>>>>>>>>>> the PnetCDF interface)
>>>>>>>>>>>>> ! cut-off domain borders
>>>>>>>>>>>>> ntx_mpi = ntx-2
>>>>>>>>>>>>> nty_mpi = nty-2
>>>>>>>>>>>>> nz_mpi  = nz
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! convert date and time integers into string 'YYYYMMDD
>>>>>>>>>>>>> hh:mm:ss'
>>>>>>>>>>>>> write (simdatetime,10)
>>>>>>>>>>>>> sim_d_YYYY,sim_d_MM,sim_d_DD,sim_t_hh,sim_t_mm,sim_t_ss
>>>>>>>>>>>>> 10 format (I4,'-',I2,'-',I2,1X,I2,':',I2,':',I2)
>>>>>>>>>>>>> simdate = simdatetime(1:10)
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! switch between cartesian axes/time (default) and  
>>>>>>>>>>>>> geographic
>>>>>>>>>>>>> axes/date (linear approx)
>>>>>>>>>>>>> coord_scale    = 1._real4
>>>>>>>>>>>>> coord_x_offset = 0._real4
>>>>>>>>>>>>> coord_x_name   = 'x  '
>>>>>>>>>>>>> coord_x_lname  = 'x_coor   '
>>>>>>>>>>>>> coord_x_units  = 'm            '
>>>>>>>>>>>>> coord_y_offset = 0._real4
>>>>>>>>>>>>> coord_y_name   = 'y  '
>>>>>>>>>>>>> coord_y_lname  = 'y_coor   '
>>>>>>>>>>>>> coord_y_units  = 'm            '
>>>>>>>>>>>>> coord_z_name   = 'z  '
>>>>>>>>>>>>> coord_t_offset = 0._real4
>>>>>>>>>>>>> coord_t_name   = 'time'
>>>>>>>>>>>>> coord_t_units  = 'sec                              '
>>>>>>>>>>>>> if (grid2geo) then
>>>>>>>>>>>>> call
>>>>>>>>>>>>> gridconvert 
>>>>>>>>>>>>> (coord_scale 
>>>>>>>>>>>>> ,coord_x_offset,coord_x_name,coord_x_lname,coord_x_units,
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> &
>>>>>>>>>>>>>
>>>>>>>>>>>>> coord_y_offset,coord_y_name,coord_y_lname,coord_y_units, &
>>>>>>>>>>>>>
>>>>>>>>>>>>> coord_z_name,coord_t_offset,coord_t_name,coord_t_units)
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! === create netCDF dataset: enter define mode
>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>> OP = 'enter define mode'
>>>>>>>>>>>>> call handle_err(nfmpi_create(my_cart, trim(outfile), 0,
>>>>>>>>>>>>> MPI_INFO_NULL,NCID),OP)
>>>>>>>>>>>>> ncid_def = NCID
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! === create netCDF dataset: assign metadata
>>>>>>>>>>>>> ----------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! ---                        set global attributes
>>>>>>>>>>>>> ----------------!
>>>>>>>>>>>>> OP = 'set global attributes'
>>>>>>>>>>>>> if (out_type == 1) then
>>>>>>>>>>>>> toffset = timetot - spinup
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,  
>>>>>>>>>>>>> 'title',
>>>>>>>>>>>>> 13_MPI_OFFSET_KIND, 'atham_picture'),OP)
>>>>>>>>>>>>> else
>>>>>>>>>>>>> toffset = timetot - spinup - nrep*periodt
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,  
>>>>>>>>>>>>> 'title',
>>>>>>>>>>>>> 11_MPI_OFFSET_KIND, 'atham_movie'),OP)
>>>>>>>>>>>>> end if
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,  
>>>>>>>>>>>>> 'content',
>>>>>>>>>>>>> 223_MPI_OFFSET_KIND, &
>>>>>>>>>>>>> 'This data corresponds to a simulation run of the
>>>>>>>>>>>>> mesoscale-beta/gamma type non-hydrostatic Active Tracer
>>>>>>>>>>>>> High-resolution Atmospheric Model ATHAM ' &
>>>>>>>>>>>>> // '(Oberhuber et al., 1998, Herzog, 1998), for a setup as
>>>>>>>>>>>>> specified under history'),OP) ! ######### add further
>>>>>>>>>>>>> description as
>>>>>>>>>>>>> necessary
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,  
>>>>>>>>>>>>> 'version',
>>>>>>>>>>>>> 21_MPI_OFFSET_KIND, 'no current versioning'),OP) !  
>>>>>>>>>>>>> #########
>>>>>>>>>>>>> change
>>>>>>>>>>>>> once available
>>>>>>>>>>>>>
>>>>>>>>>>>>> call date_and_time(date, time)
>>>>>>>>>>>>> call getlog(login)
>>>>>>>>>>>>> if (volcano_setup) then
>>>>>>>>>>>>> isetup=1
>>>>>>>>>>>>> else if (coignimbrite_setup) then
>>>>>>>>>>>>> isetup=2
>>>>>>>>>>>>> else if (nuclear_setup) then
>>>>>>>>>>>>> isetup=3
>>>>>>>>>>>>> else if (convection_setup) then
>>>>>>>>>>>>> isetup=4
>>>>>>>>>>>>> else if (procsconfig_setup) then
>>>>>>>>>>>>> isetup=5
>>>>>>>>>>>>> else
>>>>>>>>>>>>> isetup=6
>>>>>>>>>>>>> end if
>>>>>>>>>>>>> hist_string = trim(date//' '//time//'  '// 
>>>>>>>>>>>>> trim(login)//'  atham
>>>>>>>>>>>>> '//model_setup(isetup))
>>>>>>>>>>>>>
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>> 'history',
>>>>>>>>>>>>> 63_MPI_OFFSET_KIND, hist_string),OP) ! 1 line/ 
>>>>>>>>>>>>> modification:  -
>>>>>>>>>>>>> cmd
>>>>>>>>>>>>> arg
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>> 'Conventions', 28_MPI_OFFSET_KIND, 'COARDS,
>>>>>>>>>>>>> //cf-pcmdi.llnl.gov/'),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_text(NCID, NF_GLOBAL,
>>>>>>>>>>>>> 'sim_onset_date_time', 19_MPI_OFFSET_KIND,  
>>>>>>>>>>>>> simdatetime),OP)
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_real(NCID, NF_GLOBAL,
>>>>>>>>>>>>> 'sim_time_offset', NF_FLOAT, 1_MPI_OFFSET_KIND,
>>>>>>>>>>>>> real(toffset,real4)),OP)
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_real(NCID, NF_GLOBAL,
>>>>>>>>>>>>> 'x_zoom_location', NF_FLOAT, 1_MPI_OFFSET_KIND,
>>>>>>>>>>>>> real(x(max(1,iindex)),real4)),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! ---                        define dimensions
>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>> OP = 'define dimensions'
>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, trim(coord_x_name),
>>>>>>>>>>>>> ntx_mpi,             x_dimID),OP)
>>>>>>>>>>>>> if (nty>4) then
>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, trim(coord_y_name),
>>>>>>>>>>>>> nty_mpi,             y_dimID),OP)
>>>>>>>>>>>>> else
>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, trim(coord_y_name),
>>>>>>>>>>>>> 1_MPI_OFFSET_KIND,   y_dimID),OP)
>>>>>>>>>>>>> end if
>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, trim(coord_z_name),
>>>>>>>>>>>>> nz_mpi,              z_dimID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_dim(NCID, coord_t_name,
>>>>>>>>>>>>> nfmpi_unlimited,     t_dimID),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! ---                        define variables
>>>>>>>>>>>>> ---------------------!
>>>>>>>>>>>>> ! dimensions
>>>>>>>>>>>>> OP = 'define variables - dimensions'
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, trim(coord_x_name),
>>>>>>>>>>>>> NF_FLOAT, 1, x_dimID, x_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, trim(coord_y_name),
>>>>>>>>>>>>> NF_FLOAT, 1, y_dimID, y_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, trim(coord_z_name),
>>>>>>>>>>>>> NF_FLOAT, 1, z_dimID, z_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, coord_t_name,
>>>>>>>>>>>>> NF_FLOAT, 1, t_dimID, t_varID),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! dynamics data
>>>>>>>>>>>>> OP = 'define variables - data'
>>>>>>>>>>>>> if (out_type == 1 .OR. DynMovie_flag == .true.) then
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'u',       NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), u_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'v',       NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), v_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'w',       NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), w_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'pnew',    NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), pnew_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'tetnew',  NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), tetnew_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'tempnew', NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), tempnew_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'density', NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), density_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'turbhor', NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), turbhor_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'turbver', NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), turbver_varID),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_def_var(NCID, 'turblen', NF_FLOAT,  
>>>>>>>>>>>>> 4,
>>>>>>>>>>>>> (/x_dimID,y_dimID,z_dimID,t_dimID/), turblen_varID),OP)
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! tracer data
>>>>>>>>>>>>> if (out_type == 1) then
>>>>>>>>>>>>> do i=1,ntgas_picture
>>>>>>>>>>>>>  call handle_err(nfmpi_def_var(NCID, var_tgas_picture(i),
>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>> tgas_varID(i)),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>  call handle_err(nfmpi_def_var(NCID, var_trac_picture(i),
>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>> trac_varID(i)),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntpas_picture
>>>>>>>>>>>>>  call handle_err(nfmpi_def_var(NCID, var_tpas_picture(i),
>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>> tpas_varID(i)),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>  var='d_'//trim(var_trac_picture(i))
>>>>>>>>>>>>>  call handle_err(nfmpi_def_var(NCID, var,
>>>>>>>>>>>>> NF_FLOAT, 3, (/x_dimID,y_dimID,t_dimID/),  
>>>>>>>>>>>>> d_trac_varID(i)),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> else if (out_type == 2) then
>>>>>>>>>>>>> do i=1,ntgas_movie
>>>>>>>>>>>>>  call handle_err(nfmpi_def_var(NCID, var_tgas_movie(i),
>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>> tgas_varID(i)),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntrac_movie
>>>>>>>>>>>>>  call handle_err(nfmpi_def_var(NCID, var_trac_movie(i),
>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>> trac_varID(i)),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntpas_movie
>>>>>>>>>>>>>  call handle_err(nfmpi_def_var(NCID, var_tpas_movie(i),
>>>>>>>>>>>>> NF_FLOAT, 4, (/x_dimID,y_dimID,z_dimID,t_dimID/),
>>>>>>>>>>>>> tpas_varID(i)),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> endif
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! ---                        assign attributes
>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>> ! common attributes: long_name, units, valid_range,
>>>>>>>>>>>>> scale_factor,
>>>>>>>>>>>>> add_offset, _FillValue, history, conventions
>>>>>>>>>>>>> ! dimensions
>>>>>>>>>>>>>
>>>>>>>>>>>>> OP = 'assign attributes - dimensions'
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, x_varID,
>>>>>>>>>>>>> 'long_name',               9_MPI_OFFSET_KIND,
>>>>>>>>>>>>> trim(coord_x_lname)),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, x_varID,
>>>>>>>>>>>>> 'units',                  13_MPI_OFFSET_KIND,
>>>>>>>>>>>>> trim(coord_x_units)),OP)
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_real(NCID, x_varID,
>>>>>>>>>>>>> 'scale_factor', NF_FLOAT,  1_MPI_OFFSET_KIND,  
>>>>>>>>>>>>> coord_scale),OP)
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_real(NCID, x_varID,
>>>>>>>>>>>>> 'add_offset',   NF_FLOAT,  1_MPI_OFFSET_KIND,
>>>>>>>>>>>>> coord_x_offset),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, y_varID,
>>>>>>>>>>>>> 'long_name',               9_MPI_OFFSET_KIND,
>>>>>>>>>>>>> trim(coord_y_lname)),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, y_varID,
>>>>>>>>>>>>> 'units',                  13_MPI_OFFSET_KIND,
>>>>>>>>>>>>> trim(coord_y_units)),OP)
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_real(NCID, y_varID,
>>>>>>>>>>>>> 'scale_factor', NF_FLOAT,  1_MPI_OFFSET_KIND,  
>>>>>>>>>>>>> coord_scale),OP)
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_real(NCID, y_varID,
>>>>>>>>>>>>> 'add_offset',   NF_FLOAT,  1_MPI_OFFSET_KIND,
>>>>>>>>>>>>> coord_y_offset),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, z_varID,
>>>>>>>>>>>>> 'long_name',               6_MPI_OFFSET_KIND,  
>>>>>>>>>>>>> 'height'),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, z_varID,
>>>>>>>>>>>>> 'units',                   1_MPI_OFFSET_KIND, 'm'),OP)
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_text(NCID, z_varID,
>>>>>>>>>>>>> 'positive',                2_MPI_OFFSET_KIND, 'up'),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, t_varID,
>>>>>>>>>>>>> 'long_name',               4_MPI_OFFSET_KIND,  
>>>>>>>>>>>>> coord_t_name),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_text(NCID, t_varID,
>>>>>>>>>>>>> 'units',                  33_MPI_OFFSET_KIND,
>>>>>>>>>>>>> trim(coord_t_units)),OP)
>>>>>>>>>>>>> !    call handle_err(nfmpi_put_att_real(NCID, t_varID,
>>>>>>>>>>>>> 'add_offset',   NF_FLOAT,  1_MPI_OFFSET_KIND,
>>>>>>>>>>>>> coord_t_offset),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !    OP = 'exit define mode NOW for testing
>>>>>>>>>>>>> ****************************************************' !
>>>>>>>>>>>>> ah_testing
>>>>>>>>>>>>> ********************************
>>>>>>>>>>>>> !    call
>>>>>>>>>>>>> handle_err(nfmpi_enddef(NCID),OP)
>>>>>>>>>>>>> ! ah_testing ********************************
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! dynamics data
>>>>>>>>>>>>> OP = 'assign attributes - data'
>>>>>>>>>>>>> if (out_type == 1 .OR. DynMovie_flag == .true.) then
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID, u_varID,
>>>>>>>>>>>>> 'long_name',              6_MPI_OFFSET_KIND, 'x_wind'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID, u_varID,
>>>>>>>>>>>>> 'units',                  7_MPI_OFFSET_KIND, 'm  
>>>>>>>>>>>>> sec-1'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID, u_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID, u_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>> call handle_err(nfmpi_put_att_real(NCID, u_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID, v_varID,
>>>>>>>>>>>>> 'long_name',              6_MPI_OFFSET_KIND, 'y_wind'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID, v_varID,
>>>>>>>>>>>>> 'units',                  7_MPI_OFFSET_KIND, 'm  
>>>>>>>>>>>>> sec-1'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID, v_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID, v_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID, v_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID, w_varID,
>>>>>>>>>>>>> 'long_name',             19_MPI_OFFSET_KIND,
>>>>>>>>>>>>> 'upward_air_velocity'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID, w_varID,
>>>>>>>>>>>>> 'units',                  7_MPI_OFFSET_KIND, 'm  
>>>>>>>>>>>>> sec-1'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID, w_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID, w_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID, w_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>> 'long_name',             20_MPI_OFFSET_KIND,
>>>>>>>>>>>>> 'air_pressure_anomaly'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>> 'units',                  2_MPI_OFFSET_KIND, 'Pa'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> pnew_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>> 'long_name',             25_MPI_OFFSET_KIND,
>>>>>>>>>>>>> 'air_potential_temperature'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>> 'units',                  1_MPI_OFFSET_KIND, 'K'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> tetnew_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>> 'long_name',             15_MPI_OFFSET_KIND,
>>>>>>>>>>>>> 'air_temperature'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>> 'units',                  7_MPI_OFFSET_KIND,  
>>>>>>>>>>>>> 'Celsius'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, -273.15),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> tempnew_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>> 'long_name',             11_MPI_OFFSET_KIND,  
>>>>>>>>>>>>> 'air_density'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>> 'units',                  6_MPI_OFFSET_KIND, 'kg m-3'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> density_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>> 'long_name',             15_MPI_OFFSET_KIND,
>>>>>>>>>>>>> 'hor_turb_energy'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>> 'units',                  8_MPI_OFFSET_KIND, 'm2  
>>>>>>>>>>>>> sec-2'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> turbhor_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>> 'long_name',             15_MPI_OFFSET_KIND,
>>>>>>>>>>>>> 'ver_turb_energy'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>> 'units',                  8_MPI_OFFSET_KIND, 'm2  
>>>>>>>>>>>>> sec-2'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> turbver_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>> 'long_name',             17_MPI_OFFSET_KIND,
>>>>>>>>>>>>> 'turb_length_scale'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_text(NCID,  
>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>> 'units',                  1_MPI_OFFSET_KIND, 'm'),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1.),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>> !       call handle_err(nfmpi_put_att_real(NCID,  
>>>>>>>>>>>>> turblen_varID,
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! tracer data
>>>>>>>>>>>>> if (out_type == 1) then
>>>>>>>>>>>>> print *, 'NetCDF output type (1-picture,2-movie)  
>>>>>>>>>>>>> ',out_type !
>>>>>>>>>>>>> #####  DO NOT REMOVE!!! optimization problem, trace source
>>>>>>>>>>>>> do i=1,ntgas_picture
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, tgas_varID(i),
>>>>>>>>>>>>> 'long_name',             25_MPI_OFFSET_KIND,
>>>>>>>>>>>>> des_tgas_picture(i)),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, tgas_varID(i),
>>>>>>>>>>>>> 'units',                  9_MPI_OFFSET_KIND, 'gram  
>>>>>>>>>>>>> kg-1'),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, trac_varID(i),
>>>>>>>>>>>>> 'long_name',             25_MPI_OFFSET_KIND,
>>>>>>>>>>>>> des_trac_picture(i)),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, trac_varID(i),
>>>>>>>>>>>>> 'units',                  9_MPI_OFFSET_KIND, 'gram  
>>>>>>>>>>>>> kg-1'),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntpas_picture
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, tpas_varID(i),
>>>>>>>>>>>>> 'long_name',             25_MPI_OFFSET_KIND,
>>>>>>>>>>>>> des_tpas_picture(i)),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, tpas_varID(i),
>>>>>>>>>>>>> 'units',                  1_MPI_OFFSET_KIND, '%'),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 100.),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntrac_picture
>>>>>>>>>>>>>  des='deposited_'//trim(des_trac_picture(i))
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, d_trac_varID(i),
>>>>>>>>>>>>> 'long_name',             35_MPI_OFFSET_KIND, des),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, d_trac_varID(i),
>>>>>>>>>>>>> 'units',                  9_MPI_OFFSET_KIND, 'gram  
>>>>>>>>>>>>> m-2'),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, d_trac_varID(i),
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, d_trac_varID(i),
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, d_trac_varID(i),
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> else if (out_type == 2) then
>>>>>>>>>>>>> print *, 'NetCDF output type (1-picture,2-movie)  
>>>>>>>>>>>>> ',out_type !
>>>>>>>>>>>>> #####  DO NOT REMOVE!!! optimization problem, trace source
>>>>>>>>>>>>> do i=1,ntgas_movie
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, tgas_varID(i),
>>>>>>>>>>>>> 'long_name',             25_MPI_OFFSET_KIND,
>>>>>>>>>>>>> des_tgas_movie(i)),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, tgas_varID(i),
>>>>>>>>>>>>> 'units',                  9_MPI_OFFSET_KIND, 'gram  
>>>>>>>>>>>>> kg-1'),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tgas_varID(i),
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntrac_movie
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, trac_varID(i),
>>>>>>>>>>>>> 'long_name',             25_MPI_OFFSET_KIND,
>>>>>>>>>>>>> des_trac_movie(i)),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, trac_varID(i),
>>>>>>>>>>>>> 'units',                  9_MPI_OFFSET_KIND, 'gram  
>>>>>>>>>>>>> kg-1'),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 1000.),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, trac_varID(i),
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> do i=1,ntpas_movie
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, tpas_varID(i),
>>>>>>>>>>>>> 'long_name',             25_MPI_OFFSET_KIND,
>>>>>>>>>>>>> des_tpas_movie(i)),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_text(NCID, tpas_varID(i),
>>>>>>>>>>>>> 'units',                  1_MPI_OFFSET_KIND, '%'),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>> 'scale_factor',NF_FLOAT,  1_MPI_OFFSET_KIND, 100.),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>> 'add_offset',  NF_FLOAT,  1_MPI_OFFSET_KIND, 0.0),OP)
>>>>>>>>>>>>>  call handle_err(nfmpi_put_att_real(NCID, tpas_varID(i),
>>>>>>>>>>>>> '_FillValue',  NF_FLOAT,  1_MPI_OFFSET_KIND, undef),OP)
>>>>>>>>>>>>> enddo
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! === create netCDF dataset: exit define mode, enter  
>>>>>>>>>>>>> data mode
>>>>>>>>>>>>> ----!
>>>>>>>>>>>>> OP = 'exit define mode'
>>>>>>>>>>>>> call handle_err(nfmpi_enddef(NCID),OP)
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! ---                        write data values
>>>>>>>>>>>>> --------------------!
>>>>>>>>>>>>> ! dimensions (only one process writes the whole  
>>>>>>>>>>>>> dimension, i.e.
>>>>>>>>>>>>> accross all processes)
>>>>>>>>>>>>> OP = 'write data values - dimensions'
>>>>>>>>>>>>>
>>>>>>>>>>>>> call handle_err(nfmpi_begin_indep_data(NCID))
>>>>>>>>>>>>>
>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID, x_varID,
>>>>>>>>>>>>> (/1_MPI_OFFSET_KIND/), (/ 
>>>>>>>>>>>>> ntx_mpi/),real(xg(2:ntx-1),real4)),OP)
>>>>>>>>>>>>> if (nty>4) then
>>>>>>>>>>>>>  call handle_err(nfmpi_put_vara_real(NCID, y_varID,
>>>>>>>>>>>>> (/1_MPI_OFFSET_KIND/), (/nty_mpi/),
>>>>>>>>>>>>> real(yg(2:nty-1),real4)),OP)
>>>>>>>>>>>>> else
>>>>>>>>>>>>>  call handle_err(nfmpi_put_vara_real(NCID, y_varID,
>>>>>>>>>>>>> (/1_MPI_OFFSET_KIND/), (/1_MPI_OFFSET_KIND/),
>>>>>>>>>>>>> real(yg(nty-1),real4)),OP)
>>>>>>>>>>>>> end if
>>>>>>>>>>>>> call handle_err(nfmpi_put_vara_real(NCID, z_varID,
>>>>>>>>>>>>> (/1_MPI_OFFSET_KIND/), (/nz_mpi/), real(z,real4)),OP)
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> call handle_err(nfmpi_end_indep_data(NCID))
>>>>>>>>>>>>>
>>>>>>>>>>>>> contains
>>>>>>>>>>>>> !
>>>>>>>>>>>>> ----------------------------------------------------------------
>>>>>>>>>>>>>
>>>>>>>>>>>>> subroutine
>>>>>>>>>>>>> gridconvert 
>>>>>>>>>>>>> (cd_scale,cd_x_offset,cd_x_name,cd_x_lname,cd_x_units,
>>>>>>>>>>>>>
>>>>>>>>>>>>> &
>>>>>>>>>>>>>
>>>>>>>>>>>>> cd_y_offset,cd_y_name,cd_y_lname,cd_y_units, &
>>>>>>>>>>>>>
>>>>>>>>>>>>> cd_z_name,cd_t_offset,cd_t_name,cd_t_units)
>>>>>>>>>>>>> ! IF: cd_* = coordinate parameters
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! linear conversion from cartesian grid into lon/lat
>>>>>>>>>>>>> system    !
>>>>>>>>>>>>> ! only valid as first approxiamtion over small
>>>>>>>>>>>>> domains         !
>>>>>>>>>>>>> ! to be replaced once a dyn core in a lat/lon system
>>>>>>>>>>>>> available !
>>>>>>>>>>>>>
>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only             : deglat_start,
>>>>>>>>>>>>> deglon_start
>>>>>>>>>>>>>
>>>>>>>>>>>>> use phys_constants, only           : pi
>>>>>>>>>>>>>
>>>>>>>>>>>>> real(real4),        intent(inout) :: cd_scale,  
>>>>>>>>>>>>> cd_x_offset,
>>>>>>>>>>>>> cd_y_offset, cd_t_offset
>>>>>>>>>>>>> character(len=3),   intent(inout) :: cd_x_name, cd_y_name,
>>>>>>>>>>>>> cd_z_name
>>>>>>>>>>>>> character(len=4),   intent(inout) :: cd_t_name
>>>>>>>>>>>>> character(len=9),   intent(inout) :: cd_x_lname,  
>>>>>>>>>>>>> cd_y_lname
>>>>>>>>>>>>> character(len=13),  intent(inout) :: cd_x_units,  
>>>>>>>>>>>>> cd_y_units
>>>>>>>>>>>>> character(len=33),  intent(inout) :: cd_t_units
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! local
>>>>>>>>>>>>> variables                                              !
>>>>>>>>>>>>>
>>>>>>>>>>>>> !--------------------------------------------------------------!
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> real(kreal)                       :: r
>>>>>>>>>>>>> real(kreal), parameter            :: a =  
>>>>>>>>>>>>> 6378137._kreal     !
>>>>>>>>>>>>> (WGS 84)
>>>>>>>>>>>>> real(kreal), parameter            :: b =  
>>>>>>>>>>>>> 6356752.3142_kreal !
>>>>>>>>>>>>> (WGS 84)
>>>>>>>>>>>>>
>>>>>>>>>>>>> ! ellipsoid radius at ATHAM domain start
>>>>>>>>>>>>> r = a*b / sqrt((b*b*(cos(deglat_start*pi/180._kreal))**2
>>>>>>>>>>>>> )+(a*a*(sin(deglat_start*pi/180._kreal))**2 ))
>>>>>>>>>>>>>
>>>>>>>>>>>>> cd_scale    = real((1/r) * (180._kreal/pi),real4)
>>>>>>>>>>>>>
>>>>>>>>>>>>> cd_x_offset = deglon_start
>>>>>>>>>>>>> cd_x_name   = 'lon'
>>>>>>>>>>>>> cd_x_lname  = 'longitude'
>>>>>>>>>>>>> cd_x_units  = 'degrees_east '
>>>>>>>>>>>>>
>>>>>>>>>>>>> cd_y_offset = deglat_start
>>>>>>>>>>>>> cd_y_name   = 'lat'
>>>>>>>>>>>>> cd_y_lname  = 'latitude '
>>>>>>>>>>>>> cd_y_units  = 'degrees_north'
>>>>>>>>>>>>>
>>>>>>>>>>>>> cd_z_name   = 'lev'
>>>>>>>>>>>>>
>>>>>>>>>>>>> coord_t_offset = real((sim_t_hh*3600 + sim_t_mm*60
>>>>>>>>>>>>> +sim_t_ss),real4)
>>>>>>>>>>>>> coord_t_name   = 'date'
>>>>>>>>>>>>> coord_t_units  = 'seconds since '//simdate//' 00:00:00
>>>>>>>>>>>>> [local]'
>>>>>>>>>>>>> ! COARDS conventions require time zone from UTC to be
>>>>>>>>>>>>> specified
>>>>>>>>>>>>> ! not necessary for local area version of ATHAM
>>>>>>>>>>>>>
>>>>>>>>>>>>> end subroutine gridconvert
>>>>>>>>>>>>> !
>>>>>>>>>>>>> ----------------------------------------------------------------
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> end subroutine netcdf_define
>>>>>>>>>>>>> ! 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> ==========================================================
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> !
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> ==========================================================
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> subroutine handle_err(status,op)
>>>>>>>>>>>>>
>>>>>>>>>>>>> use atham_module, only : myrank
>>>>>>>>>>>>>
>>>>>>>>>>>>> integer,                    intent (in   ) :: status
>>>>>>>>>>>>> character(len=*), optional, intent (in   ) :: op
>>>>>>>>>>>>>
>>>>>>>>>>>>> if ( status /= nf_noerr ) then
>>>>>>>>>>>>> if (myrank==0) then
>>>>>>>>>>>>> print *, ""
>>>>>>>>>>>>> print *, "P-NetCDF operation: ",trim(op)
>>>>>>>>>>>>> print *, trim(nfmpi_strerror(status))
>>>>>>>>>>>>> print *, "atham_pnetcdf error: program execution stopped"
>>>>>>>>>>>>> print *, ""
>>>>>>>>>>>>>
>>>>>>>>>>>>> stop
>>>>>>>>>>>>> end if
>>>>>>>>>>>>> end if
>>>>>>>>>>>>>
>>>>>>>>>>>>> return
>>>>>>>>>>>>>
>>>>>>>>>>>>> end subroutine handle_err
>>>>>>>>>>>>> !
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> = 
>>>>>>>>>>>>> ==========================================================
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> end module atham_pnetcdf
>>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>> -- 
>>>>>>>>> ______________________________________________________________________
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> Alex HOFFMANN                     PhD Candidate
>>>>>>>>> Centre for Atmospheric Science
>>>>>>>>> Department of Geography           University of Cambridge
>>>>>>>>> Downing Place, Cambridge CB2 3EN, UK
>>>>>>>>> e-mail: ah519 at cam.ac.uk  tel: +44 (0)1223 766581
>>>>>>>>> www.geog.cam.ac.uk/people/hoffmann/
>>>>>>>>> ______________________________________________________________________
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>> -- 
>>>>>>> ______________________________________________________________________
>>>>>>>
>>>>>>>
>>>>>>> Alex HOFFMANN                     PhD Candidate
>>>>>>> Centre for Atmospheric Science
>>>>>>> Department of Geography           University of Cambridge
>>>>>>> Downing Place, Cambridge CB2 3EN, UK
>>>>>>> e-mail: ah519 at cam.ac.uk  tel: +44 (0)1223 766581
>>>>>>> www.geog.cam.ac.uk/people/hoffmann/
>>>>>>> ______________________________________________________________________
>>>>>>>
>>>>>>> !-*- F90 -*-
>>>>>>> program mpi_setup_netCDFout
>>>>>>>
>>>>>>> !--------------------------------------------------------------------!
>>>>>>>
>>>>>>> ! author:  Alex
>>>>>>> HOFFMANN                                             !
>>>>>>> ! email:
>>>>>>> ah519 at cam.ac.uk                                           !
>>>>>>> ! date:    February
>>>>>>> 2009                                             !
>>>>>>> ! version:
>>>>>>> development                                               !
>>>>>>> !
>>>>>>> !
>>>>>>> ! simple program to test mpi runs & netCDF output (serial &
>>>>>>> parallel)!
>>>>>>> !--------------------------------------------------------------------!
>>>>>>>
>>>>>>>
>>>>>>> !  test code designed to be used with the variables from the
>>>>>>> mesoscale
>>>>>>> !     atmospheric model ATHAM
>>>>>>> !  re-structured as stand-alone
>>>>>>> !  reads and writes as netCDF an index matrix for a 4x4 matrix  
>>>>>>> in
>>>>>>> real4
>>>>>>>
>>>>>>> !  use precision
>>>>>>> !  use atham_module, only : nx,ny,ntx,nty,npx,npy
>>>>>>>
>>>>>>> implicit none
>>>>>>>
>>>>>>> !private
>>>>>>> !public
>>>>>>>
>>>>>>> ! --- mpi interface for Fortran 77 only, use header  
>>>>>>> files-----------
>>>>>>> include 'mpif.h'
>>>>>>>
>>>>>>>
>>>>>>> ! --- precision  
>>>>>>> ----------------------------------------------------
>>>>>>> integer, parameter         :: real4 =
>>>>>>> selected_real_kind(6),            &
>>>>>>>                           real8 =
>>>>>>> selected_real_kind(12),           &
>>>>>>>                           kreal = real8
>>>>>>> integer, parameter         :: int4 =
>>>>>>> selected_int_kind(7),              &
>>>>>>>                           int8 =
>>>>>>> selected_int_kind(13),             &
>>>>>>>                           kint = int4
>>>>>>> ! --- atham_module  
>>>>>>> -------------------------------------------------
>>>>>>> integer(kint)                            ::  
>>>>>>> nx,ny,ntx,nty,npx,npy
>>>>>>> integer(kint)                            :: my_ranks,
>>>>>>> nprocs,           &
>>>>>>>                                         my_cart, myrank,  
>>>>>>> mycoor(2)
>>>>>>> logical                                  :: mpi_initialized, mpi
>>>>>>> ! --- other  
>>>>>>> --------------------------------------------------------
>>>>>>> integer(kint)                            :: ierror
>>>>>>> real(real4), allocatable, dimension(:,:) :: array_g, array
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> mpi_initialized=.false.
>>>>>>> ! ---for initial  
>>>>>>> testing--------------------------------------------
>>>>>>> nx  = 4_kint ! user-specified number of grid points along x
>>>>>>> ny  = 4_kint
>>>>>>> npx = 2_kint ! number of processes along x
>>>>>>> npy = 2_kint
>>>>>>>
>>>>>>> ntx = nx     ! re-adjusted total number of grid points along x
>>>>>>> nty = ny
>>>>>>> nx  = nx/npx ! number of grid points along x per process
>>>>>>> ny  = ny/npy
>>>>>>>
>>>>>>>
>>>>>>> ! --- MAIN
>>>>>>> ! --- mpi setup calls  
>>>>>>> ----------------------------------------------
>>>>>>> call initialize_mpi
>>>>>>> call setup_mpi
>>>>>>>
>>>>>>> ! --- netCDF out calls  
>>>>>>> ---------------------------------------------
>>>>>>> call nc_out_s
>>>>>>> call nc_out_p
>>>>>>>
>>>>>>> ! --- mpi finalize call  
>>>>>>> --------------------------------------------
>>>>>>> call mpi_Finalize(ierror)
>>>>>>> ! --- END MAIN
>>>>>>>
>>>>>>> contains
>>>>>>> ! 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> ================================================================
>>>>>>>
>>>>>>> subroutine initialize_mpi()
>>>>>>>
>>>>>>> !    use atham_module, only : mpi, mpi_initialized
>>>>>>> !    use atham_module, only : myrank, nprocs
>>>>>>>
>>>>>>> integer(kint)              :: ierror, errorclass
>>>>>>> character(len=100)         :: message
>>>>>>>
>>>>>>> ! assign variables for serial version
>>>>>>> ierror=0
>>>>>>> myrank=0
>>>>>>> nprocs=1
>>>>>>> mpi=.false.
>>>>>>> ! initialize MPI
>>>>>>> mpi=.true.
>>>>>>> if (.not. mpi_initialized) call mpi_init(ierror)
>>>>>>> call MPI_Comm_rank(mpi_comm_world,myrank,ierror)
>>>>>>> call MPI_Comm_size(mpi_comm_world,nprocs,ierror)
>>>>>>> mpi=.true.
>>>>>>> mpi_initialized=.true.
>>>>>>>
>>>>>>> if (myrank==0) then
>>>>>>>  print *, 'Myrank 0: ',myrank
>>>>>>>  print *, 'MPI initialized!'
>>>>>>>  print *, 'Number of processes [nprocs]: ',nprocs
>>>>>>> endif
>>>>>>>
>>>>>>> if (ierror>0) then
>>>>>>>  call MPI_Error_string(ierror,message,100)
>>>>>>>  call MPI_Error_class(ierror,errorclass)
>>>>>>>  call exec_stop(message,errorclass)
>>>>>>> endif
>>>>>>>
>>>>>>> end subroutine initialize_mpi
>>>>>>> ! 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> ================================================================
>>>>>>>
>>>>>>> subroutine setup_mpi()
>>>>>>>
>>>>>>> !    use atham_module, only : mpi, mpi_initialized
>>>>>>> !    use atham_module, only : myrank,nprocs,my_cart,mycoor
>>>>>>>
>>>>>>> integer(kint)              :: ierror, errorclass
>>>>>>> character(len=100)         :: message
>>>>>>> integer(kint)              :: idims(2),shape_array(2)
>>>>>>> logical                    :: reorder
>>>>>>> integer(kint)              :: i,j,istart,jstart
>>>>>>>
>>>>>>> if (nprocs/=npx*npy) then
>>>>>>>  print *, 'nprocs /= npx*npy - execution aborted'
>>>>>>>  call mpi_abort(mpi_comm_world,ierror)
>>>>>>>  stop
>>>>>>> endif
>>>>>>> ! create virtual topology
>>>>>>> idims(1:2) = (/npx,npy/)
>>>>>>> reorder = .true.
>>>>>>>
>>>>>>> call
>>>>>>> MPI_Cart_create(mpi_comm_world,2,idims, 
>>>>>>> (/.false.,.false./),reorder,my_cart,ierror)
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> call MPI_Comm_rank(my_cart,myrank,ierror) ! ranks on new  
>>>>>>> communicator
>>>>>>> call MPI_Cart_coords(my_cart,myrank,2,mycoor,ierror)
>>>>>>>
>>>>>>> if(myrank==0) then
>>>>>>>  print *, 'Coords of process ',myrank,' :',mycoor
>>>>>>> endif
>>>>>>>
>>>>>>> if (ierror>0) then
>>>>>>>  call MPI_Error_string(ierror,message,100)
>>>>>>>  call MPI_Error_class(ierror,errorclass)
>>>>>>>  call exec_stop(message,errorclass)
>>>>>>> endif
>>>>>>>
>>>>>>> ! allocate and assign global grid on each process
>>>>>>> allocate(array_g(ntx,nty))
>>>>>>>
>>>>>>> array_g(1,:) = (/ 11.0,12.0,13.0,14.0 /)
>>>>>>> array_g(2,:) = (/ 21.0,22.0,23.0,24.0 /)
>>>>>>> array_g(3,:) = (/ 31.0,32.0,33.0,34.0 /)
>>>>>>> array_g(4,:) = (/ 41.0,42.0,43.0,44.0 /)
>>>>>>> if (myrank==0) then
>>>>>>>  print *, 'whole 4*4 array on each process:'
>>>>>>>  print *, array_g(1,:)
>>>>>>>  print *, array_g(2,:)
>>>>>>>  print *, array_g(3,:)
>>>>>>>  print *, array_g(4,:)
>>>>>>> endif
>>>>>>>
>>>>>>> ! allocate and assign subgrid portions on various processes
>>>>>>> allocate(array(nx,ny))
>>>>>>> istart = mycoor(1)*nx
>>>>>>> jstart = mycoor(2)*ny
>>>>>>>
>>>>>>> do i=1,nx
>>>>>>>  do j=1,ny
>>>>>>>     array(i,j) = array_g(istart+i,jstart+j)
>>>>>>>  enddo
>>>>>>> enddo
>>>>>>>
>>>>>>> ! output results from one process
>>>>>>> if (myrank==0) then
>>>>>>>  print *, ''
>>>>>>>  print *, 'subgrid array on process with coords ',mycoor,':'
>>>>>>>  shape_array = shape(array)
>>>>>>>  print *, 'shape: ',shape_array
>>>>>>>  do i=1,shape_array(1)
>>>>>>>     print *, array(i,:)
>>>>>>>  enddo
>>>>>>>  print *, ''
>>>>>>> endif
>>>>>>>
>>>>>>> end subroutine setup_mpi
>>>>>>> ! 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> ================================================================
>>>>>>>
>>>>>>> subroutine nc_out_s()
>>>>>>>
>>>>>>> use netcdf
>>>>>>>
>>>>>>> !    use atham_module, only : makename
>>>>>>> !    use atham_module, only : mycoor
>>>>>>> !    use atham_module, only : nx,ny
>>>>>>>
>>>>>>> integer(kint)        :: NCID
>>>>>>> character(len=10)    :: name
>>>>>>> integer(kint)        ::  
>>>>>>> x_dimID,y_dimID,x_varID,y_varID,array_varID
>>>>>>>
>>>>>>>
>>>>>>> call makename('out',mycoor,name)
>>>>>>> name = trim(name)//'.nc'
>>>>>>>
>>>>>>> call handle_err_nc_s(nf90_create(path = './'//trim(name),  
>>>>>>> cmode =
>>>>>>> 0, ncid = NCID))
>>>>>>> call handle_err_nc_s(nf90_def_dim(ncid = NCID, name = 'x',      
>>>>>>> len
>>>>>>> = nx, dimid = x_dimID))
>>>>>>> call handle_err_nc_s(nf90_def_dim(ncid = NCID, name = 'y',      
>>>>>>> len
>>>>>>> = ny, dimid = y_dimID))
>>>>>>> call handle_err_nc_s(nf90_def_var(ncid = NCID, name = 'x',
>>>>>>> xtype = NF90_FLOAT, dimids = x_dimID, varid = x_varID))
>>>>>>> call handle_err_nc_s(nf90_def_var(ncid = NCID, name = 'y',
>>>>>>> xtype = NF90_FLOAT, dimids = y_dimID, varid = y_varID))
>>>>>>> call handle_err_nc_s(nf90_def_var(ncid = NCID, name = 'array',
>>>>>>> xtype = NF90_FLOAT, dimids = (/x_dimID,y_dimID/), varid =
>>>>>>> array_varID))
>>>>>>> call handle_err_nc_s(nf90_enddef(ncid = NCID))
>>>>>>> call handle_err_nc_s(nf90_put_var(ncid = NCID, varid =  
>>>>>>> array_varID,
>>>>>>> values = array))
>>>>>>> call handle_err_nc_s(nf90_close(ncid = NCID))
>>>>>>>
>>>>>>> end subroutine nc_out_s
>>>>>>> ! 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> ================================================================
>>>>>>>
>>>>>>> subroutine nc_out_p
>>>>>>>
>>>>>>> !    use atham_module, only : mycoor,my_cart,myrank
>>>>>>> !    use atham_module, only : nx,ny
>>>>>>>
>>>>>>> ! precompiler problem: include file had to be hardwired into  
>>>>>>> code
>>>>>>> !    include 'pnetcdf.inc'
>>>>>>> #include "/usr/local/parallel-netcdf-1.0.1/include/pnetcdf.inc"
>>>>>>>
>>>>>>>
>>>>>>> integer(kint)                :: NCIDP
>>>>>>> character(len=10)            :: name
>>>>>>> integer(kint)                ::
>>>>>>> x_dimIDP,y_dimIDP,x_varIDP,y_varIDP,array_varIDP
>>>>>>> integer(kint)                :: dimIDPs(2)
>>>>>>> integer(MPI_OFFSET_KIND)     :: dimlenx,dimleny
>>>>>>> integer(MPI_OFFSET_KIND)     :: istart,jstart
>>>>>>> integer(MPI_OFFSET_KIND)     :: ntx_mpi,nty_mpi,nx_mpi,ny_mpi
>>>>>>>
>>>>>>> ntx_mpi = ntx
>>>>>>> nty_mpi = nty
>>>>>>> nx_mpi  = nx
>>>>>>> ny_mpi  = ny
>>>>>>>
>>>>>>> name   = 'out_glb.nc'
>>>>>>> istart = mycoor(1)*nx_mpi + 1
>>>>>>> jstart = mycoor(2)*ny_mpi + 1
>>>>>>>
>>>>>>>
>>>>>>> call
>>>>>>> handle_err_nc_p(nfmpi_create(my_cart,'./'//trim(name), 
>>>>>>> 0,MPI_INFO_NULL,NCIDP),0)
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> call handle_err_nc_p(nfmpi_def_dim(NCIDP,'x',ntx_mpi,x_dimIDP), 
>>>>>>> 1)
>>>>>>>
>>>>>>> call handle_err_nc_p(nfmpi_def_dim(NCIDP,'y',nty_mpi,y_dimIDP), 
>>>>>>> 2)
>>>>>>>
>>>>>>> call
>>>>>>> handle_err_nc_p(nfmpi_def_var(NCIDP,'x',NF_FLOAT, 
>>>>>>> 1,x_dimIDP,x_varIDP),3)
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> call
>>>>>>> handle_err_nc_p(nfmpi_def_var(NCIDP,'y',NF_FLOAT, 
>>>>>>> 1,y_dimIDP,y_varIDP),4)
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> call
>>>>>>> handle_err_nc_p(nfmpi_def_var(NCIDP,'array',NF_FLOAT,2,(/ 
>>>>>>> x_dimIDP,y_dimIDP/),array_varIDP),5)
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> call handle_err_nc_p(nfmpi_enddef(NCIDP),6)
>>>>>>>
>>>>>>> if (myrank==0) then
>>>>>>>  print *, ''
>>>>>>>  print *, 'array in nc_out_p subroutine'
>>>>>>>  print *, array
>>>>>>>  print *, ''
>>>>>>>  print *, '(/nx,ny/)         ',(/nx,ny/)
>>>>>>>  print *, '(/istart,jstart/) ',(/istart,jstart/)
>>>>>>> endif
>>>>>>>
>>>>>>> if (myrank==0) then
>>>>>>>  print *, ''
>>>>>>>  call handle_err_nc_p(nfmpi_inq_dimlen(NCIDP,x_dimIDP,dimlenx), 
>>>>>>> 7)
>>>>>>>  call handle_err_nc_p(nfmpi_inq_dimlen(NCIDP,y_dimIDP,dimleny), 
>>>>>>> 8)
>>>>>>>  print *, 'pnetCDF var dimensions in xy ',dimlenx,dimleny
>>>>>>>  print *, 'ntx_mpi ',ntx_mpi
>>>>>>>  print *, 'ntx     ',ntx
>>>>>>>  print *, 'nx_mpi  ',nx_mpi
>>>>>>>  print *, 'nx      ',nx
>>>>>>> endif
>>>>>>>
>>>>>>> ! this should have been possible in collective data mode, but  
>>>>>>> did
>>>>>>> not seem to work !!!
>>>>>>> call handle_err_nc_p(nfmpi_begin_indep_data(NCIDP),9)
>>>>>>>
>>>>>>> call
>>>>>>> handle_err_nc_p(nfmpi_put_vara_real(NCIDP,array_varIDP,(/ 
>>>>>>> istart,jstart/),(/nx_mpi,ny_mpi/),reshape(array,(/ 
>>>>>>> nx_mpi,ny_mpi/))),10)
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> call handle_err_nc_p(nfmpi_end_indep_data(NCIDP),11)
>>>>>>>
>>>>>>> call handle_err_nc_p(nfmpi_close(NCIDP),12)
>>>>>>>
>>>>>>> end subroutine nc_out_p
>>>>>>> ! 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> ================================================================
>>>>>>>
>>>>>>> subroutine handle_err_nc_s(status)
>>>>>>>
>>>>>>> use netcdf
>>>>>>>
>>>>>>> integer,                    intent (in   ) :: status
>>>>>>>
>>>>>>> if ( status /= nf90_noerr ) then
>>>>>>>  print *, ""
>>>>>>>  print *, trim(nf90_strerror(status))
>>>>>>>  print *, "atham_netcdf error: program execution stopped"
>>>>>>>  print *, ""
>>>>>>>  stop
>>>>>>> end if
>>>>>>> return
>>>>>>> end subroutine handle_err_nc_s
>>>>>>> ! 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> ================================================================
>>>>>>>
>>>>>>> subroutine handle_err_nc_p(status,callID)
>>>>>>>
>>>>>>> !    use atham_module, only : myrank
>>>>>>>
>>>>>>> #include "/usr/local/parallel-netcdf-1.0.1/include/pnetcdf.inc"
>>>>>>>
>>>>>>> integer,                    intent (in   ) :: status
>>>>>>> integer,                    intent (in   ) :: callID
>>>>>>>
>>>>>>> if ( status /= nf_noerr ) then
>>>>>>>  if (myrank==0) then
>>>>>>>     print *, ""
>>>>>>>     print *, trim(nfmpi_strerror(status))
>>>>>>>     print *, "atham_pnetcdf error: program execution stopped
>>>>>>> ",callID
>>>>>>>     print *, ""
>>>>>>> !         stop
>>>>>>>  end if
>>>>>>> endif
>>>>>>> return
>>>>>>> end subroutine handle_err_nc_p
>>>>>>> ! 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> ================================================================
>>>>>>>
>>>>>>> subroutine exec_stop(MESSAGE,ERRORCLASS)
>>>>>>>
>>>>>>> character(len=*), intent(in   ) :: MESSAGE
>>>>>>> integer(kint),    intent(in   ) :: ERRORCLASS
>>>>>>>
>>>>>>> integer(kint)                   :: ierror
>>>>>>>
>>>>>>> print *, 'ERROR: MPI execution stopped'
>>>>>>> print *, 'error class: ', ERRORCLASS
>>>>>>> print *, 'message: ', trim(MESSAGE)
>>>>>>> call mpi_abort(mpi_comm_world,ierror)
>>>>>>> stop
>>>>>>>
>>>>>>> end subroutine exec_stop
>>>>>>> ! 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> ================================================================
>>>>>>>
>>>>>>> subroutine makename(prefix,icoor,name)
>>>>>>> !------------------------------------------------------------------!
>>>>>>> ! create character  
>>>>>>> name=prefix_(icoor(1)+1)_(icoor(2)+1)           !
>>>>>>> !------------------------------------------------------------------!
>>>>>>> character(len=*), intent(in)             :: prefix
>>>>>>> integer(kint), dimension(2), intent(in)  :: icoor
>>>>>>> character(len=*), intent(out)            :: name
>>>>>>> !------------------------------------------------------------------!
>>>>>>> ! local  
>>>>>>> variables                                                  !
>>>>>>> !------------------------------------------------------------------!
>>>>>>> character(len=25)                        :: numx,numy
>>>>>>> integer(kint)                            :: lenx,leny,lenf
>>>>>>>
>>>>>>> write(numx,'(i5)') icoor(1)+1
>>>>>>> numx=adjustl(numx)
>>>>>>> lenx=len_trim(numx)
>>>>>>>
>>>>>>> write(numy,'(i5)') icoor(2)+1
>>>>>>> numy=adjustl(numy)
>>>>>>> leny=len_trim(numy)
>>>>>>>
>>>>>>> name=adjustl(prefix)
>>>>>>> lenf=len_trim(name)
>>>>>>>
>>>>>>> name=name(1:lenf)//'_'//numx(1:lenx)//'_'//numy(1:leny)
>>>>>>>
>>>>>>> end subroutine makename
>>>>>>> ! 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> = 
>>>>>>> ================================================================
>>>>>>>
>>>>>>>
>>>>>>> end program mpi_setup_netCDFout
>>>>>>> # MAKEFILE used for development module testing
>>>>>>> # Alex Hoffmann
>>>>>>>
>>>>>>> # set netCDF path
>>>>>>> NETCDF = /usr/local/netcdf-3.6.2-ifort
>>>>>>> PNETCDF= /usr/local/parallel-netcdf-1.0.1
>>>>>>>
>>>>>>> # set compiler
>>>>>>> F90    = mpif90  # mpif90 refers to a shell script that includes
>>>>>>> links
>>>>>>> to the libraries
>>>>>>>
>>>>>>>
>>>>>>> # set compiler flags
>>>>>>> OFLAGS = -g
>>>>>>> # set include flags
>>>>>>> IFLAGS = -I$(NETCDF)/include/ -I$(PNETCDF)/include/
>>>>>>> # set pre-processing flags
>>>>>>> #PFLAGS = -fpp
>>>>>>> # concatenate flags
>>>>>>> FLAGS  = $(OFLAGS) $(IFLAGS) $(PFLAGS)
>>>>>>> # set library flags
>>>>>>> LFLAGS = -L$(NETCDF)/lib -lnetcdf -L$(PNETCDF)/lib -lpnetcdf
>>>>>>>
>>>>>>> RM = rm -f
>>>>>>>
>>>>>>> #MAIN =
>>>>>>> MAIN = mpi_setup_netCDFout
>>>>>>> OBJECTS =
>>>>>>>
>>>>>>>
>>>>>>> all: $(MAIN)
>>>>>>>
>>>>>>>
>>>>>>> $(MAIN): $(MAIN).F90 $(OBJECTS)
>>>>>>> $(F90) $(FLAGS) -o $@ $^ $(LFLAGS)
>>>>>>>
>>>>>>> clean:
>>>>>>> $(RM) *.o *.mod
>>>>>>
>>>>>
>>>>> -- 
>>>>> ______________________________________________________________________
>>>>>
>>>>> Alex HOFFMANN                     PhD Candidate
>>>>> Centre for Atmospheric Science
>>>>> Department of Geography           University of Cambridge
>>>>> Downing Place, Cambridge CB2 3EN, UK
>>>>> e-mail: ah519 at cam.ac.uk  tel: +44 (0)1223 766581
>>>>> www.geog.cam.ac.uk/people/hoffmann/
>>>>> ______________________________________________________________________
>>>>>
>>>>
>>>
>>> -- 
>>> ______________________________________________________________________
>>>
>>> Alex HOFFMANN                     PhD Candidate
>>> Centre for Atmospheric Science
>>> Department of Geography           University of Cambridge
>>> Downing Place, Cambridge CB2 3EN, UK
>>> e-mail: ah519 at cam.ac.uk  tel: +44 (0)1223 766581
>>> www.geog.cam.ac.uk/people/hoffmann/
>>> ______________________________________________________________________
>>>
>>
>
> -- 
> ______________________________________________________________________
>
> Alex HOFFMANN                     PhD Candidate
> Centre for Atmospheric Science
> Department of Geography           University of Cambridge
> Downing Place, Cambridge CB2 3EN, UK
> e-mail: ah519 at cam.ac.uk  tel: +44 (0)1223 766581
> www.geog.cam.ac.uk/people/hoffmann/
> ______________________________________________________________________
>



More information about the parallel-netcdf mailing list