crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-27 05:35:49 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Exposing C++ create_voxel, fixing indexing bug, and removing fortran implementation.

Browse source
This commit is contained in:
Patrick Shriwise 2018-10-09 22:10:38 -05:00
parent a83ea14a86
commit 3fab6fdb2e
3 changed files with 21 additions and 145 deletions
Download patch file Download diff file Expand all files Collapse all files

2
include/openmc/plot.h
View file

@ -78,6 +78,8 @@ void output_ppm(ObjectPlot* pl,
extern "C" void create_ppm(ObjectPlot* pl);
extern "C" void create_voxel(ObjectPlot *pl);
extern "C" void position_rgb(Particle* p, ObjectPlot* pl, int rgb[3], int &id);
extern "C" void voxel_init(hid_t file_id, const hsize_t* dims, hid_t* dspace,

131
src/plot.F90
View file

@ -40,6 +40,11 @@ module plot
type(ObjectPlot), intent(in) :: pl
end subroutine create_ppm
subroutine create_voxel(pl) bind(C)
import ObjectPlot
type(ObjectPlot), intent(in) :: pl
end subroutine create_voxel
end interface
contains
@ -71,130 +76,4 @@ contains
err = 0
end function openmc_plot_geometry
!===============================================================================
! CREATE_VOXEL outputs a binary file that can be input into silomesh for 3D
! geometry visualization. It works the same way as create_ppm by dragging a
! particle across the geometry for the specified number of voxels. The first
! 3 int(4)'s in the binary are the number of x, y, and z voxels. The next 3
! real(8)'s are the widths of the voxels in the x, y, and z directions. The next
! 3 real(8)'s are the x, y, and z coordinates of the lower left point. Finally
! the binary is filled with entries of four int(4)'s each. Each 'row' in the
! binary contains four int(4)'s: 3 for x,y,z position and 1 for cell or material
! id. For 1 million voxels this produces a file of approximately 15MB.
!===============================================================================
subroutine create_voxel(pl)
type(ObjectPlot), intent(in) :: pl
integer(C_INT) :: x, y, z ! voxel location indices
integer :: rgb(3) ! colors (red, green, blue) from 0-255
integer :: id ! id of cell or material
integer(C_INT), target :: data(pl%pixels(3),pl%pixels(2))
integer(HID_T) :: file_id
integer(HID_T) :: dspace
integer(HID_T) :: memspace
integer(HID_T) :: dset
integer(HSIZE_T) :: dims(3)
real(8) :: vox(3) ! x, y, and z voxel widths
real(8) :: ll(3) ! lower left starting point for each sweep direction
type(Particle) :: p
type(ProgressBar) :: progress
interface
subroutine voxel_init(file_id, dims, dspace, dset, memspace) bind(C)
import HID_T, HSIZE_T
integer(HID_T), value :: file_id
integer(HSIZE_T), intent(in) :: dims(*)
integer(HID_T), intent(out) :: dspace
integer(HID_T), intent(out) :: dset
integer(HID_T), intent(out) :: memspace
end subroutine voxel_init
subroutine voxel_write_slice(x, dspace, dset, memspace, buf) bind(C)
import C_INT, HID_T, C_PTR
integer(C_INT), value :: x
integer(HID_T), value :: dspace
integer(HID_T), value :: dset
integer(HID_T), value :: memspace
type(C_PTR), value :: buf
end subroutine voxel_write_slice
subroutine voxel_finalize(dspace, dset, memspace) bind(C)
import HID_T
integer(HID_T), value :: dspace
integer(HID_T), value :: dset
integer(HID_T), value :: memspace
end subroutine voxel_finalize
end interface
! compute voxel widths in each direction
vox = pl % width/dble(pl % pixels)
! initial particle position
ll = pl % origin - pl % width / TWO
! allocate and initialize particle
call particle_initialize(p)
p % coord(1) % xyz = ll
p % coord(1) % uvw = [ HALF, HALF, HALF ]
p % coord(1) % universe = root_universe
! Open binary plot file for writing
file_id = file_open(pl%path_plot, 'w')
! write header info
call write_attribute(file_id, "filetype", 'voxel')
call write_attribute(file_id, "version", VERSION_VOXEL)
call write_attribute(file_id, "openmc_version", VERSION)
#ifdef GIT_SHA1
call write_attribute(file_id, "git_sha1", GIT_SHA1)
#endif
! Write current date and time
call write_attribute(file_id, "date_and_time", time_stamp())
call write_attribute(file_id, "num_voxels", pl%pixels)
call write_attribute(file_id, "voxel_width", vox)
call write_attribute(file_id, "lower_left", ll)
! Create dataset for voxel data -- note that the dimensions are reversed
! since we want the order in the file to be z, y, x
dims(:) = pl % pixels
call voxel_init(file_id, dims, dspace, dset, memspace)
! move to center of voxels
ll = ll + vox / TWO
do x = 1, pl % pixels(1)
call progress % set_value(dble(x)/dble(pl % pixels(1))*100)
do y = 1, pl % pixels(2)
do z = 1, pl % pixels(3)
! get voxel color
call position_rgb(p, pl, rgb, id)
! write to plot file
data(z,y) = id
! advance particle in z direction
p % coord(1) % xyz(3) = p % coord(1) % xyz(3) + vox(3)
end do
! advance particle in y direction
p % coord(1) % xyz(2) = p % coord(1) % xyz(2) + vox(2)
p % coord(1) % xyz(3) = ll(3)
end do
! advance particle in y direction
p % coord(1) % xyz(1) = p % coord(1) % xyz(1) + vox(1)
p % coord(1) % xyz(2) = ll(2)
p % coord(1) % xyz(3) = ll(3)
! Write to HDF5 dataset
call voxel_write_slice(x, dspace, dset, memspace, c_loc(data))
end do
call voxel_finalize(dspace, dset, memspace)
call file_close(file_id)
end subroutine create_voxel
end module plot

33
src/plot.cpp
View file

@ -376,33 +376,27 @@ void create_voxel(ObjectPlot *pl) {
#endif
// Write current date and time
// write_attribute(file_id, "date_and_time", time_stamp());
// write_attribute(file_id, "date_and_time", time_stamp()); <-- RE-ADD!!!
hsize_t three = 3;
write_attr_int(file_id, 1, &three, "num_voxels", pl->pixels);
write_attr_double(file_id, 1, &three, "voxel_width", vox);
write_attr_double(file_id, 1, &three, "lower_left", ll);
// Create dataset for voxel data -- note that the dimensions are reversed
// since we want the order in the file to be z, y, x
hsize_t dims[3];
dims[0] = pl->pixels[0];
dims[1] = pl->pixels[1];
dims[2] = pl->pixels[2];
hid_t dspace, dset, memspace;
voxel_init(file_id, &(dims[0]), &dspace, &dset, &memspace);
// move to center of voxels
ll[0] = ll[0] + vox[0] / TWO;
ll[1] = ll[1] + vox[1] / TWO;
ll[2] = ll[2] + vox[2] / TWO;
ImageData data;
data.resize(pl->pixels[2]);
for (auto & i : data) {
i.resize(pl->pixels[1]);
for (auto & j : i) {
j.resize(1);
}
}
int data[pl->pixels[1]][pl->pixels[2]];
int rgb[3], id;
for (int x = 0; x < pl->pixels[0]; x++) {
@ -411,6 +405,8 @@ void create_voxel(ObjectPlot *pl) {
for(int z = 0; z < pl->pixels[2]; z++) {
// get voxel color
position_rgb(p, pl, rgb, id);
// write to plot data
data[y][z] = id;
// advance particle in z direction
p->coord[0].xyz[2] = p->coord[0].xyz[2] + vox[2];
}
@ -419,12 +415,11 @@ void create_voxel(ObjectPlot *pl) {
p->coord[0].xyz[2] = ll[2];
}
// advance particle in x direction
p->coord[0].xyz[0] = p->coord[0].xyz[0] + vox[0];
p->coord[0].xyz[1] = ll[1];
p->coord[0].xyz[2] = ll[2];
// Write to HDF5 dataset
voxel_write_slice(x, dspace, dset, memspace, &(data[0]));
p->coord[0].xyz[0] = p->coord[0].xyz[0] + vox[0];
p->coord[0].xyz[1] = ll[1];
p->coord[0].xyz[2] = ll[2];
// Write to HDF5 dataset
voxel_write_slice(x, dspace, dset, memspace, &(data[0]));
}
voxel_finalize(dspace, dset, memspace);
@ -456,7 +451,7 @@ voxel_init(hid_t file_id, const hsize_t* dims, hid_t* dspace, hid_t* dset,
void
voxel_write_slice(int x, hid_t dspace, hid_t dset, hid_t memspace, void* buf)
{
hssize_t offset[3] {x - 1, 0, 0};
hssize_t offset[3] {x, 0, 0};
H5Soffset_simple(dspace, offset);
H5Dwrite(dset, H5T_NATIVE_INT, memspace, dspace, H5P_DEFAULT, buf);
}