Creating NEMO output diagnostics

Hi there,

I would like to create my own output diagnostics for NEMO, which could be for example, density, or Froude and Rossby numbers. But first I am trying to figure out how nemo calls and outputs current diagnostics, and I am a little confused. Hopefully someone can help.

My understanding is that the ocean diagnostics are defined in the field_def_nemo-oce.xml file. Some are then calculated or called in src/OCE/DIA/diawri.F90. For example, the square of the buoyancy frequency is just iom_put( "bn2", rn2 ) on line 392, and this works. But I cannot find where rn2 comes into diawri from. How does this work? I know rn2 is first defined in OCE/oce.F90 along with many other terms. I have tried adding in in-situ density rhd by defining it in field_def_nemo-oce.xml, adding it to file_def_nemo.xml, and adding the line iom_put( "rhd", rhd ) to diawri.F90. But this did not work.

I’m going to guess for the square of buoyancy frequency, there is a CALL bn( ... ) from the src/OCE/TRA/eosbn2.F90 to calculate the term, and somewhere along the line it gets passed in to diawri.F90?

Thanks in advance.


I think it’s a lot more simple than that, rn2 is defined in the oce module as a public variable, and passed along to diawri, that just bins it into the output files.
And for your “rhd” variable, are you sure the iom_put line is called indeed and/or that it contains
data when the call is made ?

Hope this helps,



Sorry, my ‘it does not work’ is pretty useless. The field for rhd is populated by NaN’s, so it is creating netcdf output, which is good. When rn2 is called in diawri.F90, there is no other mention of rn2 except the one on line 392. Where is the field of rn2 being populated? It must be calculated somewhere. The density rhd is seemingly not populated when I call iom_put for rhd. Both in-situ density and buoyancy frequency squared are calculated in eosbn2.F90.

I guess I want to avoid calculating terms more than I have to, since an easy solution would be to calculate rhd in diawri.F90.

Apologies if I’m missing something obvious.



Yes, I think you miss something obvious. diawri.F90 is not made to compute anything but to write things down. rn2 is calculated, there is a call to bn2 in step.F90, which fills it with data. And since it is a public variable of the oce module then diawri.F90 has access to it and just writes it down.
Regarding rhd, I have no idea of why it prints only NaN, perhaps it simply has no value when entering diawri.F90, just print like the total of rhd in an ascii file to check if it contains values or just NaNs.



Hi Robinson,

I think you’ve cleared things up for me. I now have a populated density field.


1 Like

Hi Thomas and Robinson,

I agree with what Robinson said. Two other things that might be worth mentioning, just in case:

  1. New diagnostic variables can be written from anywhere in the code, not just diawri.F90. Just use call iom_put("<variable_name_in_xml_file>", <data_array>), where <data_array> can be a shared module variable, but doesn’t have to (it can also be a temporary array). If your new diagnostic implies specific computations which are not needed for the actual run, even better is to nest it all (calculations and iom_put) into a IF iom_use conditional (or combination thereof), so that in the case where you’re not asking for that specific variable, these specific computations will not be made uselessly. The iom routines are part of the iom module, so you might need to import it if it wasn’t already.

  2. My personal experience is that it’s best not to make NEMO output any NaN’s - fill the data you’re outputting with zeros or something else at initialization. XIOS doesn’t like NaN’s (numerical conversions errors, etc.), so a crash might occur if you change ARCH file, HPC environment, etc.