Skip to content

Commit

Permalink
Merge pull request #449 from joakim-hove/grid-unit-system
Browse files Browse the repository at this point in the history
Read gridunit from grid file
  • Loading branch information
joakim-hove authored Aug 17, 2018
2 parents 2369ea7 + ec86fc9 commit afa3a91
Show file tree
Hide file tree
Showing 3 changed files with 123 additions and 16 deletions.
1 change: 1 addition & 0 deletions lib/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -347,6 +347,7 @@ foreach (name ecl_alloc_cpgrid
ecl_grid_copy
ecl_grid_create
ecl_grid_DEPTHZ
ecl_grid_unit_system
ecl_grid_export
ecl_grid_init_fwrite
ecl_grid_reset_actnum
Expand Down
75 changes: 59 additions & 16 deletions lib/ecl/ecl_grid.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -699,7 +699,7 @@ struct ecl_cell_struct {
};



static ert_ecl_unit_enum ecl_grid_check_unit_system(const ecl_kw_type * gridunit_kw);
static void ecl_grid_init_mapaxes_data_float( const ecl_grid_type * grid , float * mapaxes);
float * ecl_grid_alloc_coord_data( const ecl_grid_type * grid );
static const float * ecl_grid_get_mapaxes( const ecl_grid_type * grid );
Expand Down Expand Up @@ -1527,6 +1527,7 @@ static bool ecl_grid_alloc_cells( ecl_grid_type * grid , bool init_valid) {
*/

static ecl_grid_type * ecl_grid_alloc_empty(ecl_grid_type * global_grid,
ert_ecl_unit_enum unit_system,
int dualp_flag,
int nx,
int ny,
Expand All @@ -1553,7 +1554,7 @@ static ecl_grid_type * ecl_grid_alloc_empty(ecl_grid_type * global_grid,
grid->index_map = NULL;
grid->fracture_index_map = NULL;
grid->inv_fracture_index_map = NULL;
grid->unit_system = ECL_METRIC_UNITS;
grid->unit_system = unit_system;


if (global_grid != NULL) {
Expand Down Expand Up @@ -2416,6 +2417,7 @@ void ecl_grid_init_GRDECL_data(ecl_grid_type * ecl_grid,
*/

static ecl_grid_type * ecl_grid_alloc_GRDECL_data__(ecl_grid_type * global_grid,
ert_ecl_unit_enum unit_system,
int dualp_flag,
bool apply_mapaxes,
int nx,
Expand All @@ -2428,7 +2430,7 @@ static ecl_grid_type * ecl_grid_alloc_GRDECL_data__(ecl_grid_type * global_grid,
const int * corsnum,
int lgr_nr) {

ecl_grid_type * ecl_grid = ecl_grid_alloc_empty(global_grid , dualp_flag , nx,ny,nz,lgr_nr,true);
ecl_grid_type * ecl_grid = ecl_grid_alloc_empty(global_grid , unit_system, dualp_flag , nx,ny,nz,lgr_nr,true);
if (ecl_grid) {
if (mapaxes != NULL)
ecl_grid_init_mapaxes( ecl_grid , apply_mapaxes, mapaxes );
Expand All @@ -2452,6 +2454,7 @@ static ecl_grid_type * ecl_grid_alloc_GRDECL_data__(ecl_grid_type * global_grid,
}

static ecl_grid_type * ecl_grid_alloc_GRDECL_data__(ecl_grid_type * global_grid,
ert_ecl_unit_enum unit_system,
int dualp_flag,
bool apply_mapaxes,
int nx,
Expand All @@ -2464,7 +2467,7 @@ static ecl_grid_type * ecl_grid_alloc_GRDECL_data__(ecl_grid_type * global_grid,
const int * corsnum,
int lgr_nr) {

ecl_grid_type * ecl_grid = ecl_grid_alloc_empty(global_grid , dualp_flag , nx,ny,nz,lgr_nr,true);
ecl_grid_type * ecl_grid = ecl_grid_alloc_empty(global_grid, unit_system, dualp_flag , nx,ny,nz,lgr_nr,true);
if (ecl_grid) {
if (mapaxes != NULL)
ecl_grid_init_mapaxes( ecl_grid , apply_mapaxes, mapaxes );
Expand Down Expand Up @@ -2521,6 +2524,7 @@ static void ecl_grid_copy_content( ecl_grid_type * target_grid , const ecl_grid_

static ecl_grid_type * ecl_grid_alloc_copy__( const ecl_grid_type * src_grid, ecl_grid_type * main_grid ) {
ecl_grid_type * copy_grid = ecl_grid_alloc_empty( main_grid ,
src_grid->unit_system,
src_grid->dualp_flag ,
ecl_grid_get_nx( src_grid ) ,
ecl_grid_get_ny( src_grid ) ,
Expand Down Expand Up @@ -2625,7 +2629,10 @@ ecl_grid_type * ecl_grid_alloc_GRDECL_data(int nx,
const int * actnum,
bool apply_mapaxes,
const float * mapaxes) {

ert_ecl_unit_enum unit_system = ECL_METRIC_UNITS;
return ecl_grid_alloc_GRDECL_data__(NULL,
unit_system,
FILEHEAD_SINGLE_POROSITY,
apply_mapaxes,
nx,
Expand All @@ -2649,8 +2656,9 @@ namespace ecl {
const int * actnum,
bool apply_mapaxes,
const float * mapaxes) {

ert_ecl_unit_enum unit_system = ECL_METRIC_UNITS;
return ecl_grid_alloc_GRDECL_data__(NULL,
unit_system,
FILEHEAD_SINGLE_POROSITY,
apply_mapaxes,
nx,
Expand Down Expand Up @@ -2690,11 +2698,12 @@ static ecl_grid_type * ecl_grid_alloc_GRDECL_kw__(ecl_grid_type * global_grid ,
const ecl_kw_type * gridhead_kw ,
const ecl_kw_type * zcorn_kw ,
const ecl_kw_type * coord_kw ,
const ecl_kw_type * gridunit_kw, /* Can be NULL */
const ecl_kw_type * mapaxes_kw , /* Can be NULL */
const ecl_kw_type * corsnum_kw, /* Can be NULL */
const int * actnum_data) { /* Can be NULL */
int gtype, nx,ny,nz, lgr_nr;

int gtype, nx,ny,nz, lgr_nr;
ert_ecl_unit_enum unit_system = ECL_METRIC_UNITS;
gtype = ecl_kw_iget_int(gridhead_kw , GRIDHEAD_TYPE_INDEX);
nx = ecl_kw_iget_int(gridhead_kw , GRIDHEAD_NX_INDEX);
ny = ecl_kw_iget_int(gridhead_kw , GRIDHEAD_NY_INDEX);
Expand All @@ -2717,13 +2726,17 @@ static ecl_grid_type * ecl_grid_alloc_GRDECL_kw__(ecl_grid_type * global_grid ,
const float * mapaxes_data = NULL;
const int * corsnum_data = NULL;

if (mapaxes_kw != NULL)
if (mapaxes_kw)
mapaxes_data = ecl_grid_get_mapaxes_from_kw__(mapaxes_kw);

if (corsnum_kw != NULL)
if (corsnum_kw)
corsnum_data = ecl_kw_get_int_ptr( corsnum_kw );

if (gridunit_kw)
unit_system = ecl_grid_check_unit_system(gridunit_kw);

return ecl_grid_alloc_GRDECL_data__(global_grid ,
unit_system,
dualp_flag ,
apply_mapaxes,
nx , ny , nz ,
Expand Down Expand Up @@ -2756,13 +2769,14 @@ ecl_grid_type * ecl_grid_alloc_GRDECL_kw( int nx, int ny , int nz ,

bool apply_mapaxes = true;
ecl_kw_type * gridhead_kw = ecl_grid_alloc_gridhead_kw( nx, ny, nz, 0);

ecl_kw_type * gridunit_kw = NULL;
ecl_grid_type * ecl_grid = ecl_grid_alloc_GRDECL_kw__(NULL,
FILEHEAD_SINGLE_POROSITY,
apply_mapaxes,
gridhead_kw,
zcorn_kw,
coord_kw,
gridunit_kw,
mapaxes_kw,
NULL,
actnum_data);
Expand Down Expand Up @@ -2994,6 +3008,7 @@ static ecl_grid_type * ecl_grid_alloc_EGRID__( ecl_grid_type * main_grid , const
ecl_kw_type * coord_kw = ecl_file_iget_named_kw( ecl_file , COORD_KW , grid_nr);
ecl_kw_type * corsnum_kw = NULL;
ecl_kw_type * actnum_kw = NULL;
ecl_kw_type * gridunit_kw = NULL;
ecl_kw_type * mapaxes_kw = NULL;
int dualp_flag;
int eclipse_version;
Expand Down Expand Up @@ -3031,6 +3046,9 @@ static ecl_grid_type * ecl_grid_alloc_EGRID__( ecl_grid_type * main_grid , const

if (ecl_file_has_kw( ecl_file , CORSNUM_KW))
corsnum_kw = ecl_file_iget_named_kw( ecl_file , CORSNUM_KW , 0);

if (ecl_file_has_kw(ecl_file, GRIDUNIT_KW))
gridunit_kw = ecl_file_iget_named_kw( ecl_file, GRIDUNIT_KW, 0);
}

{
Expand All @@ -3040,8 +3058,9 @@ static ecl_grid_type * ecl_grid_alloc_EGRID__( ecl_grid_type * main_grid , const
gridhead_kw ,
zcorn_kw ,
coord_kw ,
gridunit_kw,
mapaxes_kw ,
corsnum_kw,
corsnum_kw,
actnum_data);

if (ECL_GRID_MAINGRID_LGR_NR != grid_nr) ecl_grid_set_lgr_name_EGRID(ecl_grid , ecl_file , grid_nr);
Expand Down Expand Up @@ -3105,11 +3124,11 @@ ecl_grid_type * ecl_grid_alloc_EGRID(const char * grid_file, bool apply_mapaxes)



static ecl_grid_type * ecl_grid_alloc_GRID_data__(ecl_grid_type * global_grid , int num_coords , int dualp_flag , bool apply_mapaxes, int nx, int ny , int nz , int grid_nr , int coords_size , int ** coords , float ** corners , const float * mapaxes) {
static ecl_grid_type * ecl_grid_alloc_GRID_data__(ecl_grid_type * global_grid , int num_coords , ert_ecl_unit_enum unit_system, int dualp_flag , bool apply_mapaxes, int nx, int ny , int nz , int grid_nr , int coords_size , int ** coords , float ** corners , const float * mapaxes) {
if (dualp_flag != FILEHEAD_SINGLE_POROSITY)
nz = nz / 2;
{
ecl_grid_type * grid = ecl_grid_alloc_empty( global_grid , dualp_flag , nx , ny , nz , grid_nr, false);
ecl_grid_type * grid = ecl_grid_alloc_empty( global_grid , unit_system, dualp_flag , nx , ny , nz , grid_nr, false);
if (grid) {
if (mapaxes != NULL)
ecl_grid_init_mapaxes( grid , apply_mapaxes , mapaxes);
Expand All @@ -3135,8 +3154,10 @@ static ecl_grid_type * ecl_grid_alloc_GRID_data__(ecl_grid_type * global_grid ,
*/

ecl_grid_type * ecl_grid_alloc_GRID_data(int num_coords , int nx , int ny , int nz , int coords_size , int ** coords , float ** corners , bool apply_mapaxes, const float * mapaxes) {
ert_ecl_unit_enum unit_system = ECL_METRIC_UNITS;
return ecl_grid_alloc_GRID_data__( NULL ,
num_coords ,
unit_system,
FILEHEAD_SINGLE_POROSITY , /* Does currently not support to determine dualp_flag from inspection. */
apply_mapaxes,
nx , ny , nz , 0 , coords_size , coords , corners , mapaxes);
Expand Down Expand Up @@ -3190,6 +3211,7 @@ static ecl_grid_type * ecl_grid_alloc_GRID__(ecl_grid_type * global_grid , const
int nx,ny,nz;
const float * mapaxes_data = NULL;
ecl_grid_type * grid;
ert_ecl_unit_enum unit_system = ECL_METRIC_UNITS;

// 1: Fetching header data from the DIMENS keyword.
{
Expand All @@ -3209,6 +3231,8 @@ static ecl_grid_type * ecl_grid_alloc_GRID__(ecl_grid_type * global_grid , const
}
}

if ((grid_nr == 0) && (ecl_file_has_kw( ecl_file , GRIDUNIT_KW)))
unit_system = ecl_grid_check_unit_system( ecl_file_iget_named_kw(ecl_file, GRIDUNIT_KW, 0));

/*
The number of COORDS/CORNERS blocks depends on the GRIDFILE option
Expand Down Expand Up @@ -3281,7 +3305,7 @@ static ecl_grid_type * ecl_grid_alloc_GRID__(ecl_grid_type * global_grid , const
coords_size = ecl_kw_get_size( coords_kw );
}
// Create the grid:
grid = ecl_grid_alloc_GRID_data__( global_grid , num_coords , dualp_flag , apply_mapaxes, nx , ny , nz , grid_nr , coords_size , coords , corners , mapaxes_data );
grid = ecl_grid_alloc_GRID_data__( global_grid , num_coords , unit_system, dualp_flag , apply_mapaxes, nx , ny , nz , grid_nr , coords_size , coords , corners , mapaxes_data );

free( coords );
free( corners );
Expand Down Expand Up @@ -3348,7 +3372,8 @@ ecl_grid_type * ecl_grid_alloc_GRID(const char * grid_file, bool apply_mapaxes)
which case all cells will be active.
*/
ecl_grid_type * ecl_grid_alloc_regular( int nx, int ny , int nz , const double * ivec, const double * jvec , const double * kvec , const int * actnum) {
ecl_grid_type * grid = ecl_grid_alloc_empty(NULL , FILEHEAD_SINGLE_POROSITY , nx , ny , nz , 0, true);
ert_ecl_unit_enum unit_system = ECL_METRIC_UNITS;
ecl_grid_type * grid = ecl_grid_alloc_empty(NULL , unit_system, FILEHEAD_SINGLE_POROSITY , nx , ny , nz , 0, true);
if (grid) {
const double grid_offset[3] = {0,0,0};

Expand Down Expand Up @@ -3405,7 +3430,9 @@ ecl_grid_type * ecl_grid_alloc_rectangular( int nx , int ny , int nz , double dx
all cells will be active. */
ecl_grid_type * ecl_grid_alloc_dxv_dyv_dzv( int nx, int ny , int nz , const double * dxv , const double * dyv , const double * dzv , const int * actnum)
{
ert_ecl_unit_enum unit_system = ECL_METRIC_UNITS;
ecl_grid_type* grid = ecl_grid_alloc_empty(NULL,
unit_system,
FILEHEAD_SINGLE_POROSITY,
nx, ny, nz,
/*lgr_nr=*/0, /*init_valid=*/true);
Expand Down Expand Up @@ -3448,7 +3475,9 @@ ecl_grid_type * ecl_grid_alloc_dxv_dyv_dzv( int nx, int ny , int nz , const doub

ecl_grid_type * ecl_grid_alloc_dxv_dyv_dzv_depthz( int nx, int ny , int nz , const double * dxv , const double * dyv , const double * dzv , const double * depthz , const int * actnum)
{
ert_ecl_unit_enum unit_system = ECL_METRIC_UNITS;
ecl_grid_type* grid = ecl_grid_alloc_empty(NULL,
unit_system,
FILEHEAD_SINGLE_POROSITY,
nx, ny, nz,
/*lgr_nr=*/0, /*init_valid=*/true);
Expand Down Expand Up @@ -3551,8 +3580,9 @@ ecl_grid_type * ecl_grid_alloc_dxv_dyv_dzv_depthz( int nx, int ny , int nz , con
*/

ecl_grid_type * ecl_grid_alloc_dx_dy_dz_tops( int nx, int ny , int nz , const double * dx , const double * dy , const double * dz , const double * tops , const int * actnum) {

ert_ecl_unit_enum unit_system = ECL_METRIC_UNITS;
ecl_grid_type* grid = ecl_grid_alloc_empty(NULL,
unit_system,
FILEHEAD_SINGLE_POROSITY,
nx, ny, nz,
0, true);
Expand Down Expand Up @@ -6052,6 +6082,19 @@ static ecl_kw_type * ecl_grid_alloc_gridunits_kw( ert_ecl_unit_enum output_unit
return gridunits_kw;
}

static ert_ecl_unit_enum ecl_grid_check_unit_system(const ecl_kw_type * gridunit_kw) {
const char * length_unit = ecl_kw_iget_char_ptr(gridunit_kw, 0);

if (strncmp(length_unit, "FEET", 4) == 0)
return ECL_FIELD_UNITS;

if (strncmp(length_unit, "CM", 2) == 0)
return ECL_LAB_UNITS;

return ECL_METRIC_UNITS;
}


/*****************************************************************/

static float ecl_grid_output_scaling( const ecl_grid_type * grid , ert_ecl_unit_enum output_unit) {
Expand Down
63 changes: 63 additions & 0 deletions lib/ecl/tests/ecl_grid_unit_system.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,63 @@
/*
Copyright (C) 2018 Statoil ASA, Norway.
The file 'ecl_grid_unit_system.c' is part of ERT - Ensemble based Reservoir Tool.
ERT is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
ERT is distributed in the hope that it will be useful, but WITHOUT ANY
WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE.
See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
for more details.
*/
#include <stdlib.h>
#include <stdbool.h>

#include <ert/util/test_util.hpp>
#include <ert/util/test_work_area.hpp>

#include <ert/ecl/ecl_grid.hpp>
#include <ert/ecl/ecl_kw.hpp>
#include <ert/ecl/ecl_file.hpp>


void test_EGRID(const char * filename, ert_ecl_unit_enum unit_system) {
ecl_grid_type * grid = ecl_grid_alloc_rectangular(4,4,2,1,1,1,NULL);
ecl_grid_fwrite_EGRID2( grid , filename, unit_system);
ecl_grid_free( grid );


ecl_grid_type * grid2 = ecl_grid_alloc(filename);
test_assert_int_equal( ecl_grid_get_unit_system(grid2), unit_system);
ecl_grid_free(grid2);
}



void test_GRID(const char * filename, ert_ecl_unit_enum unit_system) {
ecl_grid_type * grid = ecl_grid_alloc_rectangular(4,4,2,1,1,1,NULL);
ecl_grid_fwrite_GRID2( grid , filename, unit_system);
ecl_grid_free( grid );


ecl_grid_type * grid2 = ecl_grid_alloc(filename);
test_assert_int_equal( ecl_grid_get_unit_system(grid2), unit_system);
ecl_grid_free(grid2);
}


int main(int argc, char **argv) {
test_work_area_type * work_area = test_work_area_alloc("grid_export");

test_EGRID("METRIC.EGRID", ECL_METRIC_UNITS);
test_EGRID("FIELD.EGRID", ECL_FIELD_UNITS);
test_GRID("METRIC.GRID", ECL_METRIC_UNITS);
test_GRID("FIELD.GRID", ECL_FIELD_UNITS);

test_work_area_free(work_area);
}

0 comments on commit afa3a91

Please sign in to comment.