Skip to content

Commit

Permalink
Merge pull request #50 from ClimateGlobalChange/whannah/add_exodus_to…
Browse files Browse the repository at this point in the history 
…_scrip

add exodus to scrip converter
  • Loading branch information
@paullric
paullric authored Mar 19, 2019
2 parents 9f04d10 + 7832eb0 commit 03b5f81
Show file tree
Hide file tree
Showing 5 changed files with 312 additions and 4 deletions.
3 changes: 2 additions & 1 deletion Makefile.am
View file
Original file line number Diff line number Diff line change
Expand Up @@ -114,6 +114,7 @@ CalculateDiffNorms_SOURCES = src/CalculateDiffNorms.cpp

MeshToTxt_SOURCES = src/MeshToTxt.cpp
ShpToMesh_SOURCES = src/ShpToMesh.cpp
ConvertExodusToSCRIP_SOURCES = src/ConvertExodusToSCRIP.cpp

bin_PROGRAMS = GenerateTestData \
GenerateCSMesh GenerateRLLMesh GenerateUTMMesh GenerateICOMesh \
Expand All @@ -122,7 +123,7 @@ bin_PROGRAMS = GenerateTestData \
ApplyOfflineMap GenerateOfflineMap \
CalculateDiffNorms GenerateGLLMetaData GenerateConnectivityFile \
GenerateTransposeMap CoarsenRectilinearData \
MeshToTxt ShpToMesh
MeshToTxt ShpToMesh ConvertExodusToSCRIP


# Utility target: build but don't run tests
Expand Down
95 changes: 95 additions & 0 deletions src/ConvertExodusToSCRIP.cpp
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,95 @@
///////////////////////////////////////////////////////////////////////////////
///
/// \file ConvertExodusToSCRIP.cpp
/// \author Paul Ullrich
/// \version September 18, 2015
///
/// <remarks>
/// Copyright 2000-2014 Paul Ullrich
///
/// This file is distributed as part of the Tempest source code package.
/// Permission is granted to use, copy, modify and distribute this
/// source code and its documentation under the terms of the GNU General
/// Public License. This software is provided "as is" without express
/// or implied warranty.
/// </remarks>

#include "CommandLine.h"
#include "GridElements.h"
#include "FiniteElementTools.h"
#include "Exception.h"
#include "Announce.h"

#include "netcdfcpp.h"

#include <cmath>
#include <iostream>

///////////////////////////////////////////////////////////////////////////////

int main(int argc, char** argv) {

NcError error(NcError::silent_nonfatal);

try {
// Input filename
std::string strInputFile;

// Output scrip filename
std::string strOutputFile;

// Parse the command line
BeginCommandLine()
CommandLineString(strInputFile, "in", "");
CommandLineString(strOutputFile, "out", "");

ParseCommandLine(argc, argv);
EndCommandLine(argv)

AnnounceBanner();

// Check file names
if (strInputFile == "") {
_EXCEPTIONT("No input file specified");
}
if (strOutputFile == "") {
_EXCEPTIONT("No output file specified");
}

//---------------------------------------------------------------------------
//---------------------------------------------------------------------------

// Load input mesh
std::cout << std::endl;
std::cout << "..Loading input mesh" << std::endl;

Mesh meshIn(strInputFile);
meshIn.RemoveZeroEdges(); // Do we need this?

// Generate new mesh
// std::cout << "..Generating SCRIP format data" << std::endl;
// Mesh meshOut;

// Write the mesh
std::cout << "..Writing mesh" << std::endl;
NcFile::FileFormat eOutputFormat = NcFile::Offset64Bits;
meshIn.WriteScrip(strOutputFile, eOutputFormat);

std::cout << "..Done writing" << std::endl;


// Announce
std::cout << "..Mesh converter exited successfully" << std::endl;
std::cout << "=========================================================";
std::cout << std::endl;

return (0);

} catch(Exception & e) {
Announce(e.ToString().c_str());
return (-1);

} catch(...) {
return (-2);
}
}
207 changes: 205 additions & 2 deletions src/GridElements.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,7 @@
#include <ctime>
#include <cmath>
#include <cstring>
#include <algorithm>
#include "netcdfcpp.h"

#include "triangle.h"
Expand Down Expand Up @@ -768,6 +769,208 @@ void Mesh::Write(

///////////////////////////////////////////////////////////////////////////////

void Mesh::WriteScrip(
const std::string & strFile,
NcFile::FileFormat eFileFormat
) const {
const int ParamLenString = 33;

// Temporarily change error reporting
NcError error_temp(NcError::verbose_fatal);

//---------------------------------------------------------------------------
// Determine block sizes
std::vector<int> vecBlockSizes;
std::vector<int> vecBlockSizeFaces;
{
std::map<int, int> mapBlockSizes;
std::map<int, int>::iterator iterBlockSize;
int iBlock;
char szBuffer[ParamLenString];

for (int i = 0; i < faces.size(); i++) {
iterBlockSize = mapBlockSizes.find(faces[i].edges.size());

if (iterBlockSize == mapBlockSizes.end()) {
mapBlockSizes.insert(
std::pair<int,int>(faces[i].edges.size(), 1));
} else {
(iterBlockSize->second)++;
}
}

vecBlockSizes.resize(mapBlockSizes.size());
vecBlockSizeFaces.resize(mapBlockSizes.size());

AnnounceStartBlock("Nodes per element");
iterBlockSize = mapBlockSizes.begin();
iBlock = 1;
for (; iterBlockSize != mapBlockSizes.end(); iterBlockSize++) {
vecBlockSizes[iBlock-1] = iterBlockSize->first;
vecBlockSizeFaces[iBlock-1] = iterBlockSize->second;

Announce("Block %i (%i nodes): %i",
iBlock, vecBlockSizes[iBlock-1], vecBlockSizeFaces[iBlock-1]);

iBlock++;
}
AnnounceEndBlock(NULL);
}
//---------------------------------------------------------------------------
// Output to a NetCDF SCRIP file
NcFile ncOut(strFile.c_str(), NcFile::Replace, NULL, 0, eFileFormat);
if (!ncOut.is_valid()) {
_EXCEPTION1("Unable to open grid file \"%s\" for writing",
strFile.c_str());
}
// Find max number of corners oer all faces
int nElementCount = faces.size();
int nCornersMax = 0;
for (int i=0; i<nElementCount; i++) {
nCornersMax = std::max( nCornersMax, (int)(faces[i].edges.size()) );
}
// SCRIP dimensions
NcDim * dimGridSize = ncOut.add_dim("grid_size", nElementCount);
NcDim * dimGridCorner = ncOut.add_dim("grid_corners", nCornersMax);
NcDim * dimGridRank = ncOut.add_dim("grid_rank", 1);
// Global attributes
ncOut.add_att("api_version", 5.00f);
ncOut.add_att("version", 5.00f);
ncOut.add_att("floating_point_word_size", 8);
ncOut.add_att("file_size", 0);
//---------------------------------------------------------------------------
// Grid Area
{
NcVar * varArea = ncOut.add_var("grid_area", ncDouble, dimGridSize);
if (varArea == NULL) {
_EXCEPTIONT("Error creating variable \"grid_area\"");
}
DataArray1D<double> area(nElementCount);
for (int i=0; i<nElementCount; i++) {
area[i] = static_cast<double>( CalculateFaceArea(faces[i], nodes) );
}
varArea->set_cur((long)0);
varArea->put(area, nElementCount);
varArea->add_att("units", "radians^2");
}
//---------------------------------------------------------------------------
// Grid center and corner coordinates
{
NcVar * varCenterLat = ncOut.add_var("grid_center_lat", ncDouble, dimGridSize);
NcVar * varCenterLon = ncOut.add_var("grid_center_lon", ncDouble, dimGridSize);
NcVar * varCornerLat = ncOut.add_var("grid_corner_lat", ncDouble, dimGridSize, dimGridCorner);
NcVar * varCornerLon = ncOut.add_var("grid_corner_lon", ncDouble, dimGridSize, dimGridCorner);
if (varCenterLat == NULL) {
_EXCEPTIONT("Error creating variable \"grid_center_lat\"");
}
if (varCenterLon == NULL) {
_EXCEPTIONT("Error creating variable \"grid_center_lon\"");
}
if (varCornerLat == NULL) {
_EXCEPTIONT("Error creating variable \"grid_corner_lat\"");
}
if (varCornerLon == NULL) {
_EXCEPTIONT("Error creating variable \"grid_corner_lon\"");
}

DataArray1D<double> centerLat(nElementCount);
DataArray1D<double> centerLon(nElementCount);
DataArray2D<double> cornerLat(nElementCount, nCornersMax);
DataArray2D<double> cornerLon(nElementCount, nCornersMax);
Node corner(0,0,0);
for (int i=0; i<nElementCount; i++) {
Node center(0,0,0);
// int nCorners = faces[i].edges.size()+1;
int nCorners = faces[i].edges.size();
for (int j=0; j<nCorners; ++j) {
corner = nodes[ faces[i][j] ];
XYZtoRLL_Deg(
corner.x, corner.y, corner.z,
cornerLon[i][j],
cornerLat[i][j]);
center = center + corner;
}
center = center / nCorners;
double dMag = sqrt(center.x * center.x +
center.y * center.y +
center.z * center.z);
center.x /= dMag;
center.y /= dMag;
center.z /= dMag;
XYZtoRLL_Deg(
center.x, center.y, center.z,
centerLon[i],
centerLat[i]);
// Adjust corner logitudes
double lonDiff;
for (int j=0; j<nCorners; ++j) {
// First check for polar point
if (cornerLat[i][j]==90. || cornerLat[i][j]==-90.) {
cornerLon[i][j] = centerLon[i];
}
// Next check for corners that wrap around prime meridian
lonDiff = centerLon[i] - cornerLon[i][j];
if (lonDiff>180) {
cornerLon[i][j] = cornerLon[i][j] + (double)360.0;
}
if (lonDiff<-180) {
cornerLon[i][j] = cornerLon[i][j] - (double)360.0;
}
}
}
varCenterLat->set_cur((long)0);
varCenterLat->put(centerLat, nElementCount);
varCenterLat->add_att("units", "degrees");
varCenterLat->add_att("_FillValue", 9.96920996838687e+36 );

varCenterLon->set_cur((long)0);
varCenterLon->put(centerLon, nElementCount);
varCenterLon->add_att("units", "degrees");
varCenterLon->add_att("_FillValue", 9.96920996838687e+36 );

for (int i=0; i<nElementCount; i++) {
varCornerLat->set_cur(i,0);
varCornerLat->put(cornerLat[i], 1, nCornersMax);
varCornerLon->set_cur(i,0);
varCornerLon->put(cornerLon[i], 1, nCornersMax);
}
varCornerLat->add_att("units", "degrees");
varCornerLon->add_att("units", "degrees");
varCornerLat->add_att("_FillValue", 9.96920996838687e+36 );
varCornerLon->add_att("_FillValue", 9.96920996838687e+36 );
}
//---------------------------------------------------------------------------
// Grid mask
{
NcVar * varMask = ncOut.add_var("grid_imask", ncDouble, dimGridSize);
if (varMask == NULL) {
_EXCEPTIONT("Error creating variable \"grid_imask\"");
}
DataArray1D<double> mask(nElementCount);
for (int i = 0; i < nElementCount; i++) {
mask[i] = static_cast<double>( 1 );
}
varMask->set_cur((long)0);
varMask->put(mask, nElementCount);
varMask->add_att("_FillValue", 9.96920996838687e+36 );
}
//---------------------------------------------------------------------------
// Grid dims
{
NcVar * varDims = ncOut.add_var("grid_dims", ncInt, dimGridRank);
if (varDims == NULL) {
_EXCEPTIONT("Error creating variable \"grid_dims\"");
}
DataArray1D<int> rank(1);
rank(0) = 1;
varDims->set_cur((long)0);
varDims->put(rank, 1);
}
//---------------------------------------------------------------------------
}

///////////////////////////////////////////////////////////////////////////////

void Mesh::Read(const std::string & strFile) {

const int ParamFour = 4;
Expand Down Expand Up @@ -1171,8 +1374,8 @@ void Mesh::Read(const std::string & strFile) {
}
}

// Remove coincident nodes.
RemoveCoincidentNodes();
// Remove coincident nodes.
RemoveCoincidentNodes();
}
}

Expand Down
5 changes: 5 additions & 0 deletions src/GridElements.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -732,6 +732,11 @@ class Mesh {
NcFile::FileFormat eFileFormat = NcFile::Classic
) const;

void WriteScrip(
const std::string & strFile,
NcFile::FileFormat eFileFormat = NcFile::Classic
) const;

/// <summary>
/// Read the mesh from a NetCDF file.
/// </summary>
Expand Down
6 changes: 5 additions & 1 deletion src/Makefile.gmake
View file
Original file line number Diff line number Diff line change
Expand Up @@ -72,6 +72,8 @@ CalculateDiffNorms_FILES= CalculateDiffNorms.cpp
MeshToTxt_FILES= MeshToTxt.cpp
ShpToMesh_FILES= ShpToMesh.cpp

ConvertExodusToSCRIP_FILES= ConvertExodusToSCRIP.cpp

########################################################################
# All executables

Expand All @@ -92,7 +94,8 @@ EXEC_TARGETS= ApplyOfflineMap \
GenerateTransposeMap \
GenerateVolumetricMesh \
MeshToTxt \
ShpToMesh
ShpToMesh \
ConvertExodusToSCRIP

########################################################################
# Build rules.
Expand All @@ -119,6 +122,7 @@ CoarsenRectilinearData_EXE: $(CoarsenRectilinearData_FILES:%.cpp=$(BUILDDIR)/%.o
CalculateDiffNorms_EXE: $(CalculateDiffNorms_FILES:%.cpp=$(BUILDDIR)/%.o)
MeshToTxt_EXE: $(MeshToTxt_FILES:%.cpp=$(BUILDDIR)/%.o)
ShpToMesh_EXE: $(ShpToMesh_FILES:%.cpp=$(BUILDDIR)/%.o)
ConvertExodusToSCRIP_EXE: $(ConvertExodusToSCRIP_FILES:%.cpp=$(BUILDDIR)/%.o)

$(EXEC_TARGETS): %: $(UTIL_FILES:%.cpp=$(BUILDDIR)/%.o) %_EXE
-@$(CXX) $(LDFLAGS) -o $@ $(UTIL_FILES:%.cpp=$(BUILDDIR)/%.o) $($*_FILES:%.cpp=$(BUILDDIR)/%.o) $(LIBRARIES)
Expand Down

0 comments on commit 03b5f81

Please sign in to comment.