bdnugget/gochem
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

SDF molfile parser

Browse Source
This commit is contained in:
chemistry-software
2024-09-26 13:50:26 +02:00
parent 4b189ea78f
commit 0841f1190b
15 changed files with 1444 additions and 10 deletions
Download Patch File Download Diff File Expand all files Collapse all files

196
c/elementdata.c Normal file
View File

@ -0,0 +1,196 @@
// Periodic Table CVS file from
// https://github.com/Bowserinator/Periodic-Table-JSON
#include <ctype.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#define MAX_ELEMENTS 120
#define MAX_CSV_LINE_LENGTH 2048
typedef struct {
char name[50];
char appearance[100];
double atomic_mass;
double boil;
char category[50];
double density;
char discovered_by[50];
double melt;
double molar_heat;
char named_by[50];
int number;
int period;
int group;
char phase[10];
char source[100];
char bohr_model_image[200];
char bohr_model_3d[200];
char spectral_img[200];
char summary[500];
char symbol[3];
int xpos;
int ypos;
int wxpos;
int wypos;
char shells[50];
char electron_configuration[50];
char electron_configuration_semantic[50];
double electron_affinity;
double electronegativity_pauling;
char ionization_energies[200];
char cpk_hex[10];
char block[2];
char image_title[200];
char image_url[200];
char image_attribution[200];
} Element;
// Function to trim leading and trailing whitespaces
char *trimWhitespace(char *str) {
char *end;
// Trim leading space
while (isspace((unsigned char)*str))
str++;
if (*str == 0) // All spaces?
return str;
// Trim trailing space
end = str + strlen(str) - 1;
while (end > str && isspace((unsigned char)*end))
end--;
// Write new null terminator character
end[1] = '\0';
return str;
}
// Function to parse a CSV line considering quoted fields
void parseCSVLine(char *line, char **fields, int field_count) {
int field_index = 0;
int in_quotes = 0;
char *field_start = line;
for (char *p = line; *p; p++) {
if (*p == '\"') {
in_quotes = !in_quotes; // Toggle quote state
} else if (*p == ',' && !in_quotes) {
*p = '\0';
fields[field_index++] = trimWhitespace(field_start);
field_start = p + 1;
}
}
fields[field_index] = trimWhitespace(field_start);
}
// Function to parse the CSV file
Element *parseCSV(const char *filename, int *element_count) {
FILE *file = fopen(filename, "r");
if (!file) {
perror("Unable to open file");
return NULL;
}
Element *elements = (Element *)malloc(MAX_ELEMENTS * sizeof(Element));
if (!elements) {
perror("Unable to allocate memory");
fclose(file);
return NULL;
}
char line[MAX_CSV_LINE_LENGTH];
fgets(line, sizeof(line), file); // Skip the header line
*element_count = 0;
while (fgets(line, sizeof(line), file)) {
if (*element_count >= MAX_ELEMENTS)
break;
char *fields[35] = {0}; // Array to hold all fields in the CSV line
parseCSVLine(line, fields, 35);
Element elem = {0}; // Initialize all struct members to zero or empty
strncpy(elem.name, fields[0], sizeof(elem.name) - 1);
strncpy(elem.appearance, fields[1], sizeof(elem.appearance) - 1);
elem.atomic_mass = atof(fields[2]);
elem.boil = atof(fields[3]);
strncpy(elem.category, fields[4], sizeof(elem.category) - 1);
elem.density = atof(fields[5]);
strncpy(elem.discovered_by, fields[6], sizeof(elem.discovered_by) - 1);
elem.melt = atof(fields[7]);
elem.molar_heat = atof(fields[8]);
strncpy(elem.named_by, fields[9], sizeof(elem.named_by) - 1);
elem.number = atoi(fields[10]);
elem.period = atoi(fields[11]);
elem.group = atoi(fields[12]);
strncpy(elem.phase, fields[13], sizeof(elem.phase) - 1);
strncpy(elem.source, fields[14], sizeof(elem.source) - 1);
strncpy(elem.bohr_model_image, fields[15],
sizeof(elem.bohr_model_image) - 1);
strncpy(elem.bohr_model_3d, fields[16], sizeof(elem.bohr_model_3d) - 1);
strncpy(elem.spectral_img, fields[17], sizeof(elem.spectral_img) - 1);
strncpy(elem.summary, fields[18], sizeof(elem.summary) - 1);
strncpy(elem.symbol, fields[19], sizeof(elem.symbol) - 1);
elem.xpos = atoi(fields[20]);
elem.ypos = atoi(fields[21]);
elem.wxpos = atoi(fields[22]);
elem.wypos = atoi(fields[23]);
strncpy(elem.shells, fields[24], sizeof(elem.shells) - 1);
strncpy(elem.electron_configuration, fields[25],
sizeof(elem.electron_configuration) - 1);
strncpy(elem.electron_configuration_semantic, fields[26],
sizeof(elem.electron_configuration_semantic) - 1);
elem.electron_affinity = atof(fields[27]);
elem.electronegativity_pauling = atof(fields[28]);
strncpy(elem.ionization_energies, fields[29],
sizeof(elem.ionization_energies) - 1);
strncpy(elem.cpk_hex, fields[30], sizeof(elem.cpk_hex) - 1);
strncpy(elem.block, fields[31], sizeof(elem.block) - 1);
strncpy(elem.image_title, fields[32], sizeof(elem.image_title) - 1);
strncpy(elem.image_url, fields[33], sizeof(elem.image_url) - 1);
strncpy(elem.image_attribution, fields[34],
sizeof(elem.image_attribution) - 1);
elements[*element_count] = elem;
(*element_count)++;
}
fclose(file);
return elements;
}
int testParser() {
int element_count;
Element *elements = parseCSV("resources/periodictable.csv", &element_count);
if (elements) {
for (int i = 0; i < element_count; i++) {
printf("Element: %s (%s)\n", elements[i].name, elements[i].symbol);
printf("Atomic Mass: %.2f\n", elements[i].atomic_mass);
printf("Category: %s\n", elements[i].category);
printf("Atomic Number: %d\n", elements[i].number);
printf("Period: %d\n", elements[i].period);
printf("Group: %d\n", elements[i].group);
printf("Phase: %s\n", elements[i].phase);
printf("Summary: %s\n", elements[i].summary);
printf("Shells: %s\n", elements[i].shells);
printf("Electron Configuration: %s\n",
elements[i].electron_configuration);
printf("Electron Configuration Semantic: %s\n",
elements[i].electron_configuration_semantic);
printf("Electron Affinity: %.2f\n", elements[i].electron_affinity);
printf("Electronegativity Pauling: %.2f\n",
elements[i].electronegativity_pauling);
printf("Ionization Energies: %s\n", elements[i].ionization_energies);
printf("CPK Hex Color: %s\n", elements[i].cpk_hex);
printf("Block: %s\n", elements[i].block);
}
free(elements);
}
return 0;
}

234
c/main.c Normal file
View File

@ -0,0 +1,234 @@
#include "debug.h"
#include "elementdata.c"
#include "raylib.h"
#include "sdfparse.c"
#include <complex.h>
#include <stdbool.h>
#include <stdio.h>
#include <string.h>
//------------------------------------------------------------------------------------------
// Types and Structures Definition
//------------------------------------------------------------------------------------------
void drawAtom(Atom atom, Model sphereModel);
void drawBond(Atom atom);
static Mesh GenMeshCustom(void);
//------------------------------------------------------------------------------------
// Program main entry point
//------------------------------------------------------------------------------------
int main(int argc, char *argv[]) {
// Initialization
//--------------------------------------------------------------------------------------
const int screenWidth = 1200;
const int screenHeight = 800;
// Molecule mol = parseSDF("./resources/mol.sdf");
// Molecule mol = parseSDF("./resources/methyl-vinyl-ketone.sdf");
Molecule mol = parseSDF("./resources/sildenafil.sdf");
// normalizeCoordinates(&mol, screenWidth, screenHeight);
int element_count;
Element *elements = parseCSV("resources/periodictable.csv", &element_count);
D printf("Number of atoms: %d\n", mol.num_atoms);
for (int i = 0; i < mol.num_atoms; i++) {
D printf("Atom %d: %s (%.4f, %.4f, %.4f) - neighbours: %d\n", i + 1,
mol.atoms[i].atom_type, mol.atoms[i].x, mol.atoms[i].y,
mol.atoms[i].z, mol.atoms[i].num_neighbours);
}
InitWindow(screenWidth, screenHeight, "LOL dongs");
SetTargetFPS(60); // Set desired framerate (frames-per-second)
Camera3D camera = {0};
camera.position = (Vector3){10.0f, 10.0f, 10.0f}; // Set camera position
camera.target = (Vector3){5.0f, 5.0f, 5.0f}; // Set camera target
camera.up = (Vector3){0.0f, 1.0f, 0.0f}; // Set camera up vector
camera.fovy = 45.0f; // Set camera field of view
// camera.projection = CAMERA_ORTHOGRAPHIC; // Camera projection
// type
camera.projection = CAMERA_PERSPECTIVE; // Camera projection type
int cameraMode = CAMERA_THIRD_PERSON;
// Model sphereModel = LoadModelFromMesh(GenMeshCustom());
Model sphereModel = LoadModelFromMesh(GenMeshSphere(1.0f, 25, 25));
//--------------------------------------------------------------------------------------
// Main game loop
while (!WindowShouldClose()) // Detect window close button or ESC key
{
// Update
//----------------------------------------------------------------------------------
if (IsKeyPressed(KEY_ONE)) {
cameraMode = CAMERA_FREE;
camera.up = (Vector3){0.0f, 1.0f, 0.0f}; // Reset roll
}
if (IsKeyPressed(KEY_TWO)) {
cameraMode = CAMERA_FIRST_PERSON;
camera.up = (Vector3){0.0f, 1.0f, 0.0f}; // Reset roll
}
if (IsKeyPressed(KEY_THREE)) {
cameraMode = CAMERA_THIRD_PERSON;
camera.up = (Vector3){0.0f, 1.0f, 0.0f}; // Reset roll
}
if (IsKeyPressed(KEY_FOUR)) {
cameraMode = CAMERA_ORBITAL;
camera.up = (Vector3){0.0f, 1.0f, 0.0f}; // Reset roll
}
if (IsKeyPressed('Z'))
camera.target = (Vector3){0.0f, 0.0f, 0.0f};
UpdateCamera(&camera, cameraMode);
// UpdateCamera(&camera, CAMERA_ORBITAL);
//----------------------------------------------------------------------------------
// Draw
//----------------------------------------------------------------------------------
BeginDrawing();
ClearBackground(RAYWHITE);
BeginMode3D(camera);
// Draw atoms
for (int i = 0; i < mol.num_atoms; i++) {
drawAtom(mol.atoms[i], sphereModel);
drawBond(mol.atoms[i]);
}
EndMode3D();
// Print atom coordinates
DrawText("Atom Coordinates:", 20, 20, 20, DARKGRAY);
for (int i = 0; i < mol.num_atoms; i++) {
DrawText(TextFormat("Atom %d: %s (%.4f, %.4f, %.4f) neighbours: %d",
mol.atoms[i].idx, mol.atoms[i].atom_type,
mol.atoms[i].x, mol.atoms[i].y, mol.atoms[i].z,
mol.atoms[i].num_neighbours),
20, 40 + i * 10, 8, DARKGRAY);
// Convenient place to reset the bond drawn flags
mol.atoms[i].bonds_drawn = false;
}
EndDrawing();
//----------------------------------------------------------------------------------
}
// De-Initialization
//--------------------------------------------------------------------------------------
free(elements);
UnloadModel(sphereModel);
CloseWindow(); // Close window and OpenGL context
//--------------------------------------------------------------------------------------
return EXIT_SUCCESS;
}
void drawBond(Atom atom) {
Color bondColor = GREEN;
for (int j = 0; j < atom.num_neighbours; j++) {
if (!atom.bonds_drawn) {
if (atom.bond_orders[j] == 2) {
bondColor = BLUE;
} else if (atom.bond_orders[j] >= 3) {
bondColor = RED;
}
DrawCylinderWiresEx((Vector3){atom.x, atom.y, atom.z},
(Vector3){atom.neighbours[j]->x,
atom.neighbours[j]->y,
atom.neighbours[j]->z},
0.1f, 0.1f, 20, bondColor);
}
atom.bonds_drawn = true;
}
}
void drawAtom(Atom atom, Model sphereModel) {
Color color;
float radius;
int scalingFactor = 100;
// radius is empirical from
// https://en.wikipedia.org/wiki/Atomic_radii_of_the_elements_(data_page)
if (strcmp(atom.atom_type, "C") == 0) {
color = BLACK;
radius = 70.0f / scalingFactor;
} else if (strcmp(atom.atom_type, "O") == 0) {
color = RED;
radius = 60.0f / scalingFactor;
} else if (strcmp(atom.atom_type, "H") == 0) {
color = LIGHTGRAY;
radius = 25.0f / scalingFactor;
} else if (strcmp(atom.atom_type, "S") == 0) {
color = YELLOW;
radius = 100.0f / scalingFactor;
} else if (strcmp(atom.atom_type, "N") == 0) {
color = BLUE;
radius = 65.0f / scalingFactor;
} else {
color = GRAY;
radius = 25.0f / scalingFactor;
}
DrawModelWires(sphereModel, (Vector3){atom.x, atom.y, atom.z}, radius, color);
}
// Generate a simple triangle mesh from code
static Mesh GenMeshCustom(void)
{
Mesh mesh = { 0 };
mesh.triangleCount = 1;
mesh.vertexCount = mesh.triangleCount*3;
mesh.vertices = (float *)MemAlloc(mesh.vertexCount*3*sizeof(float)); // 3 vertices, 3 coordinates each (x, y, z)
mesh.texcoords = (float *)MemAlloc(mesh.vertexCount*2*sizeof(float)); // 3 vertices, 2 coordinates each (x, y)
mesh.normals = (float *)MemAlloc(mesh.vertexCount*3*sizeof(float)); // 3 vertices, 3 coordinates each (x, y, z)
// Vertex at (0, 0, 0)
mesh.vertices[0] = 0;
mesh.vertices[1] = 0;
mesh.vertices[2] = 0;
mesh.normals[0] = 0;
mesh.normals[1] = 1;
mesh.normals[2] = 0;
mesh.texcoords[0] = 0;
mesh.texcoords[1] = 0;
// Vertex at (1, 0, 2)
mesh.vertices[3] = 1;
mesh.vertices[4] = 0;
mesh.vertices[5] = 2;
mesh.normals[3] = 0;
mesh.normals[4] = 1;
mesh.normals[5] = 0;
mesh.texcoords[2] = 0.5f;
mesh.texcoords[3] = 1.0f;
// Vertex at (2, 0, 0)
mesh.vertices[6] = 2;
mesh.vertices[7] = 0;
mesh.vertices[8] = 0;
mesh.normals[6] = 0;
mesh.normals[7] = 1;
mesh.normals[8] = 0;
mesh.texcoords[4] = 1;
mesh.texcoords[5] =0;
// Upload mesh data from CPU (RAM) to GPU (VRAM) memory
UploadMesh(&mesh, false);
return mesh;
}

127
c/sdfparse.c Normal file
View File

@ -0,0 +1,127 @@
#include "debug.h"
#include <ctype.h>
#include <stdbool.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#define MAX_ATOMS 1000
#define MAX_LINE_LENGTH 256
#define MAX_NEIGHBOURS 4
typedef struct Atom {
float x, y, z;
char atom_type[3];
int idx;
int num_neighbours;
int bond_orders[MAX_NEIGHBOURS];
struct Atom *neighbours[MAX_NEIGHBOURS];
bool bonds_drawn;
} Atom;
typedef struct Molecule {
Atom atoms[MAX_ATOMS];
int num_atoms;
} Molecule;
void addBond(Atom *atom1, Atom *atom2, int bond_order);
Molecule parseSDF(const char *filename) {
FILE *file = fopen(filename, "r");
if (!file) {
fprintf(stderr, "Error opening file %s\n", filename);
exit(1);
}
Molecule mol;
mol.num_atoms = 0;
char line[MAX_LINE_LENGTH];
while (fgets(line, MAX_LINE_LENGTH, file) != NULL) {
// break if you reach M END
if (strncmp(line, "M END", 6) == 0) {
D printf("Reached END\n");
break;
}
// continue of line contains V2000 or V3000
if (strstr(line, "V2000") || strstr(line, "V3000"))
continue;
// Parse the atoms and their coordinates
if (isdigit(line[4]) && isdigit(line[6]) && isalpha(line[31])) {
sscanf(line, "\t%f\t%f\t%f\t%s", &mol.atoms[mol.num_atoms].x,
&mol.atoms[mol.num_atoms].y, &mol.atoms[mol.num_atoms].z,
mol.atoms[mol.num_atoms].atom_type);
mol.atoms[mol.num_atoms].idx = mol.num_atoms + 1;
mol.atoms[mol.num_atoms].bonds_drawn = false;
++mol.num_atoms;
}
// Parse the connectivity info block, atom count in SDF starts at 1
if (isdigit(line[2]) && isdigit(line[20])) {
int atom_idx, neighbour_idx, bond_order;
sscanf(line, "%d %d %d", &atom_idx, &neighbour_idx, &bond_order);
D printf("idx %d, buurman %d, multipliciteit %d\n", atom_idx,
neighbour_idx, bond_order);
addBond(&mol.atoms[atom_idx - 1], &mol.atoms[neighbour_idx - 1],
bond_order);
}
}
if (fclose(file)) {
fprintf(stderr, "Error CLOSING file %s for some reason!\n", filename);
exit(1);
};
return mol;
}
void addBond(Atom *atom1, Atom *atom2, int bond_order) {
if (atom1->num_neighbours < MAX_NEIGHBOURS &&
atom2->num_neighbours < MAX_NEIGHBOURS) {
atom1->neighbours[atom1->num_neighbours] = atom2;
atom1->bond_orders[atom1->num_neighbours] = bond_order;
atom1->num_neighbours++;
atom2->neighbours[atom2->num_neighbours] = atom1;
atom2->bond_orders[atom2->num_neighbours] = bond_order;
atom2->num_neighbours++;
} else {
D printf("ERROR: MAX_NEIGHBOURS exceeded!\n");
}
}
void normalizeCoordinates(Molecule *mol, int width, int height) {
float min_x = mol->atoms[0].x, max_x = mol->atoms[0].x;
float min_y = mol->atoms[0].y, max_y = mol->atoms[0].y;
float min_z = mol->atoms[0].z, max_z = mol->atoms[0].z;
// Find min and max coords
for (int i = 1; i < mol->num_atoms; i++) {
if (mol->atoms[i].x < min_x)
min_x = mol->atoms[i].x;
if (mol->atoms[i].x > max_x)
max_x = mol->atoms[i].x;
if (mol->atoms[i].y < min_y)
min_y = mol->atoms[i].y;
if (mol->atoms[i].y > max_y)
max_y = mol->atoms[i].y;
if (mol->atoms[i].z < min_z)
min_z = mol->atoms[i].z;
if (mol->atoms[i].z > max_z)
max_z = mol->atoms[i].z;
}
// Scaling factors
float scale_x = (float)width / (max_x - min_x) / 100;
float scale_y = (float)height / (max_y - min_y) / 100;
// Yeah well there is no depth to a screen lol
float scale_z = (float)height / (max_z - min_z) / 100;
// Normalize coordinates, lets hope we never need the original coords lmao
for (int i = 0; i < mol->num_atoms; i++) {
mol->atoms[i].x = (mol->atoms[i].x - min_x) * scale_x;
mol->atoms[i].y = (mol->atoms[i].y - min_y) * scale_y;
mol->atoms[i].z = (mol->atoms[i].z - min_z) * scale_z;
}
}