Converting SELFIES Strings to 2D Molecular Images

import selfies as sf
import random
from rdkit import Chem
from rdkit.Chem import rdMolDescriptors

# Get the complete SELFIES alphabet
alphabet = sf.get_semantic_robust_alphabet()
print(f'SELFIES alphabet contains {len(alphabet)} symbols')
print('First 20 symbols:', list(alphabet)[:20])

# Generate completely random SELFIES strings
random.seed(42)  # For reproducibility
for i in range(5):
    # Pick 8 random symbols from the alphabet
    random_selfies = ''.join(random.choices(list(alphabet), k=8))
    
    # Decode to SMILES and validate
    smiles = sf.decoder(random_selfies)
    mol = Chem.MolFromSmiles(smiles)
    formula = rdMolDescriptors.CalcMolFormula(mol) if mol else 'Invalid'
    
    print(f'{i+1}. SELFIES: {random_selfies}')
    print(f'   SMILES:  {smiles}')
    print(f'   Formula: {formula}')
    print(f'   Valid:   {mol is not None}')
    print()

import random
from rdkit import Chem

# Common SMILES characters including bonds and rings
smiles_chars = ['C', 'N', 'O', 'S', 'P', 'F', 'Cl', 'Br', 
               '(', ')', '[', ']', '=', '#', '1', '2', '3', '4', '5']

random.seed(42)  # Same seed for fair comparison
valid_count = 0

for i in range(20):
    # Generate random 8-character SMILES string
    random_smiles = ''.join(random.choices(smiles_chars, k=8))
    mol = Chem.MolFromSmiles(random_smiles)
    valid = mol is not None
    
    if valid:
        valid_count += 1
        print(f'{i+1}. SMILES: {random_smiles} ✓ VALID')
    else:
        print(f'{i+1}. SMILES: {random_smiles} ✗ Invalid')

print(f'Success rate: {valid_count}/20 ({valid_count/20*100:.1f}%)')

#!/usr/bin/env python3
import selfies as sf
from rdkit import Chem
from rdkit.Chem import Draw, rdDepictor, rdMolDescriptors
from PIL import Image, ImageDraw, ImageFont
from pathlib import Path

def selfies_to_png(selfies_string: str, output_file: str, size: int = 500):
    """Convert a SELFIES string to a PNG image with molecular formula legend."""
    
    # Input validation
    if not selfies_string or not selfies_string.strip():
        raise ValueError("SELFIES string cannot be empty")
    
    if size <= 0:
        raise ValueError(f"Image size must be positive, got: {size}")

    # Ensure output directory exists
    output_path = Path(output_file)
    output_path.parent.mkdir(parents=True, exist_ok=True)
    
    # Decode SELFIES to SMILES with comprehensive error handling
    try:
        smiles_string = sf.decoder(selfies_string.strip())
    except Exception as e:
        raise ValueError(f"Invalid SELFIES string: '{selfies_string}'. SELFIES decoding error: {e}")
    
    if not smiles_string:
        raise ValueError(f"SELFIES string '{selfies_string}' decoded to empty SMILES.")
    
    # Create RDKit molecule
    mol = Chem.MolFromSmiles(smiles_string)
    if mol is None:
        raise ValueError(f"Invalid SELFIES: decoded to invalid SMILES '{smiles_string}'")
    
    # Generate the image with legend
    img = create_molecule_image(mol, selfies_string, size)
    img.save(output_file, "PNG", optimize=True)
    
def create_molecule_image(mol: Chem.Mol, selfies_string: str, size: int = 500):
    """Create molecule image with SELFIES string in legend and dynamic sizing."""
    
    # Calculate dynamic sizes based on image size
    sizes = _calculate_dynamic_sizes(size)
    
    # Generate 2D coordinates and molecular properties
    rdDepictor.Compute2DCoords(mol)
    molecular_formula = rdMolDescriptors.CalcMolFormula(mol)
    
    # Create base molecular image
    mol_img = Draw.MolToImage(mol, size=(size, size))
    if mol_img.mode != "RGBA":
        mol_img = mol_img.convert("RGBA")
    
    # Calculate required legend height and width for long SELFIES strings
    legend_height = _calculate_legend_height(molecular_formula, selfies_string, sizes, size)
    legend_height = max(legend_height, int(size * 0.15))  # Ensure minimum height
    
    # Calculate if we need extra width for the SELFIES text
    final_width = _calculate_required_width(selfies_string, sizes, size)
    
    # Create final image with appropriate dimensions
    total_height = size + legend_height
    final_img = Image.new("RGBA", (final_width, total_height), "white")
    
    # Center the molecule image if we made the canvas wider
    mol_x_offset = (final_width - size) // 2
    final_img.paste(mol_img, (mol_x_offset, 0))
    
    # Add formatted legend with molecular formula and SELFIES
    draw = ImageDraw.Draw(final_img)
    font_regular, font_small = _load_fonts(sizes['regular_font_size'], sizes['subscript_font_size'])
    
    _draw_molecular_formula(draw, molecular_formula, font_regular, font_small, sizes, size, mol_x_offset)
    _draw_selfies_legend(draw, selfies_string, font_regular, sizes, molecular_formula, size, final_width)
    
    return final_img

def _draw_wrapped_text(
    draw: ImageDraw.Draw, text: str, start_x: int, start_y: int, max_width: int, font, sizes: dict
) -> None:
    """Draw text with word wrapping at bracket boundaries for SELFIES strings."""
    # For SELFIES, we want to break at bracket boundaries to maintain readability
    # Split SELFIES into tokens (brackets and their contents)
    tokens = []
    current_token = ""
    
    for char in text:
        if char == '[':
            if current_token:
                tokens.append(current_token)
                current_token = ""
            current_token = char
        elif char == ']':
            current_token += char
            tokens.append(current_token)
            current_token = ""
        else:
            current_token += char
    
    if current_token:
        tokens.append(current_token)
    
    # Now draw tokens, wrapping as needed at logical boundaries
    x_pos = start_x
    y_pos = start_y
    line_height = sizes['regular_font_size'] + 4  # Generous line spacing
    
    for token in tokens:
        token_width = int(draw.textlength(token, font=font))
        
        # If this token would exceed the line, wrap to next line
        if x_pos + token_width > start_x + max_width and x_pos > start_x:
            y_pos += line_height
            x_pos = start_x
        
        draw.text((x_pos, y_pos), token, fill="black", font=font)
        x_pos += token_width

def _calculate_dynamic_sizes(image_size: int):
    """Calculate dynamic sizing values based on image size for responsive design."""
    return {
        'legend_height': int(image_size * 0.08),
        'legend_y_offset': int(image_size * 0.02),
        'legend_x_offset': int(image_size * 0.02),
        'subscript_y_offset': int(image_size * 0.006),
        'regular_font_size': int(image_size * 0.028),
        'subscript_font_size': int(image_size * 0.02),
    }

# Ethanol - simplest alcohol
python selfies2png.py "[C][C][O]" -o ethanol_selfies.png --size 600

# Acetone - demonstrates branching
python selfies2png.py "[C][C][Branch1][C][=O][C]" -o acetone_selfies.png --size 600

import selfies as sf

# A molecule with branching and rings - nicotine
selfies_str = "[C][N][C][Branch1][C][C][C][C][Ring1][=Branch1][C][=C][C][N][=C][Ring1][Branch1]"
smiles_str, attribution = sf.decoder(selfies_str, attribute=True)

print(f"Original SELFIES: {selfies_str}")
print(f"Decoded SMILES:   {smiles_str}")
print("\nHow each SMILES token was built:")
for i, attr_map in enumerate(attribution):
    contributing_tokens = [a.token for a in attr_map.attribution]
    print(f"  #{i+1}: '{attr_map.token}' ← {contributing_tokens}")

import selfies as sf
from rdkit import Chem
from rdkit.Chem import Draw

# Create a hypervalent iodine compound (common in organic synthesis)
hypervalent_smiles = 'O=I(O)(O)(O)(O)O'  # Periodic acid
mol = Chem.MolFromSmiles(hypervalent_smiles)

if mol:
    # Convert to SELFIES with different constraint levels
    standard_selfies = sf.encoder(hypervalent_smiles, strict=True)
    relaxed_selfies = sf.encoder(hypervalent_smiles, strict=False)
    
    print(f"Original SMILES: {hypervalent_smiles}")
    print(f"Standard SELFIES: {standard_selfies}")
    print(f"Relaxed SELFIES:  {relaxed_selfies}")
    
    # Now decode with hypervalent constraints enabled
    sf.set_semantic_constraints("hypervalent")
    decoded_hypervalent = sf.decoder(relaxed_selfies)
    
    print(f"Hypervalent decode: {decoded_hypervalent}")
    
    # Visualize the hypervalent structure
    mol_hypervalent = Chem.MolFromSmiles(decoded_hypervalent)
    if mol_hypervalent:
        img = Draw.MolToImage(mol_hypervalent, size=(300, 300))
        img.save("hypervalent_iodine.png")

import selfies as sf
import numpy as np
from collections import Counter

# Example: preparing a small drug dataset
drug_smiles = [
    "CC(=O)OC1=CC=CC=C1C(=O)O",  # Aspirin
    "CN1C=NC2=C1C(=O)N(C(=O)N2C)C",  # Caffeine  
    "CC(C)CC1=CC=C(C=C1)C(C)C(=O)O",  # Ibuprofen
    "O=C(O)C1=CC=CC=C1",  # Benzoic acid
]

# Convert to SELFIES
drug_selfies = [sf.encoder(smiles) for smiles in drug_smiles]
print("SELFIES representations:")
for i, selfies_str in enumerate(drug_selfies):
    print(f"{i+1}. {selfies_str}")

# Build vocabulary from our dataset
alphabet = sf.get_alphabet_from_selfies(drug_selfies)
alphabet.add("[nop]")  # Padding token for variable-length sequences
alphabet = sorted(alphabet)

print(f"\nVocabulary size: {len(alphabet)}")
print(f"Tokens: {alphabet[:10]}...")  # Show first 10 tokens

# Create encoding dictionaries
token_to_idx = {token: i for i, token in enumerate(alphabet)}
idx_to_token = {i: token for token, i in token_to_idx.items()}

# Encode one molecule for neural network input
caffeine_selfies = drug_selfies[1]
tokens = list(sf.split_selfies(caffeine_selfies))
print(f"\nCaffeine tokens: {tokens}")

# Convert to numerical representation
indices = [token_to_idx[token] for token in tokens]
print(f"As indices: {indices}")

# Create padded one-hot encoding
max_len = 25  # Pad to fixed length
label_encoding, one_hot = sf.selfies_to_encoding(
    caffeine_selfies, 
    vocab_stoi=token_to_idx, 
    pad_to_len=max_len,
    enc_type="both"
)

print(f"Label encoding shape: {label_encoding.shape}")
print(f"One-hot encoding shape: {one_hot.shape}")

import selfies as sf
import random
from rdkit import Chem
from rdkit.Chem import Descriptors

def generate_random_drug_like_molecules(n=5):
    """Generate random molecules and filter for drug-like properties."""
    alphabet = list(sf.get_semantic_robust_alphabet())
    
    # Focus on common drug atoms
    drug_atoms = [token for token in alphabet if any(atom in token for atom in ['C', 'N', 'O', 'S', 'F'])]
    
    molecules = []
    for i in range(n):
        # Generate random 12-token SELFIES
        random_selfies = ''.join(random.choices(drug_atoms, k=12))
        
        # Decode and analyze
        smiles = sf.decoder(random_selfies)
        mol = Chem.MolFromSmiles(smiles)
        
        # Calculate drug-like properties
        mw = Descriptors.MolWt(mol)
        logp = Descriptors.MolLogP(mol)
        hbd = Descriptors.NumHDonors(mol)
        hba = Descriptors.NumHAcceptors(mol)
        
        # Check Lipinski's Rule of Five
        lipinski_pass = (mw <= 500 and logp <= 5 and hbd <= 5 and hba <= 10)
        
        molecules.append({
            'selfies': random_selfies,
            'smiles': smiles,
            'mw': mw,
            'logp': logp,
            'lipinski': lipinski_pass
        })
        
        print(f"Molecule {i+1}:")
        print(f"  SELFIES: {random_selfies}")
        print(f"  SMILES:  {smiles}")
        print(f"  MW: {mw:.1f}, LogP: {logp:.1f}")
        print(f"  Lipinski compliant: {lipinski_pass}")
        print()
    
    return molecules

# Generate some random drug-like molecules
random.seed(42)
drug_molecules = generate_random_drug_like_molecules(3)

import selfies as sf
from rdkit import Chem
from rdkit.Chem import Descriptors
import random

def mutate_molecule(selfies_str, mutation_rate=0.2):
    """Mutate a SELFIES string by randomly changing tokens."""
    alphabet = list(sf.get_semantic_robust_alphabet())
    tokens = list(sf.split_selfies(selfies_str))
    
    mutated_tokens = []
    for token in tokens:
        if random.random() < mutation_rate:
            # Replace with random token
            mutated_tokens.append(random.choice(alphabet))
        else:
            mutated_tokens.append(token)
    
    return ''.join(mutated_tokens)

def explore_chemical_space(starting_molecule, generations=5):
    """Perform a random walk through chemical space."""
    current_selfies = sf.encoder(starting_molecule)
    
    print(f"Starting molecule: {starting_molecule}")
    print(f"Starting SELFIES:  {current_selfies}")
    print()
    
    for gen in range(generations):
        # Mutate the molecule
        mutated_selfies = mutate_molecule(current_selfies)
        mutated_smiles = sf.decoder(mutated_selfies)
        mol = Chem.MolFromSmiles(mutated_smiles)
        
        # Calculate properties
        if mol:
            mw = Descriptors.MolWt(mol)
            logp = Descriptors.MolLogP(mol)
            
            print(f"Generation {gen+1}:")
            print(f"  SELFIES: {mutated_selfies}")
            print(f"  SMILES:  {mutated_smiles}")
            print(f"  MW: {mw:.1f}, LogP: {logp:.2f}")
            print()
            
            current_selfies = mutated_selfies

# Start with aspirin and explore
random.seed(123)
explore_chemical_space("CC(=O)OC1=CC=CC=C1C(=O)O", generations=3)

import selfies as sf
from rdkit import Chem
from rdkit.Chem import Descriptors
from concurrent.futures import ProcessPoolExecutor
import time

def process_single_selfies(selfies_str):
    """Process one SELFIES string and return properties."""
    try:
        smiles = sf.decoder(selfies_str)
        mol = Chem.MolFromSmiles(smiles)
        if mol:
            return {
                'selfies': selfies_str,
                'smiles': smiles,
                'mw': Descriptors.MolWt(mol),
                'logp': Descriptors.MolLogP(mol),
                'valid': True
            }
    except Exception as e:
        return {'selfies': selfies_str, 'valid': False, 'error': str(e)}

def benchmark_processing(selfies_list, use_parallel=True):
    """Compare serial vs parallel processing."""
    
    if use_parallel:
        # Parallel processing
        start_time = time.time()
        with ProcessPoolExecutor(max_workers=4) as executor:
            results = list(executor.map(process_single_selfies, selfies_list))
        parallel_time = time.time() - start_time
        print(f"Parallel processing: {parallel_time:.2f}s for {len(selfies_list)} molecules")
        return results
    else:
        # Serial processing  
        start_time = time.time()
        results = [process_single_selfies(s) for s in selfies_list]
        serial_time = time.time() - start_time
        print(f"Serial processing: {serial_time:.2f}s for {len(selfies_list)} molecules")
        return results

# Generate test dataset
alphabet = list(sf.get_semantic_robust_alphabet())
test_selfies = [''.join(random.choices(alphabet, k=10)) for _ in range(100)]

# Benchmark both approaches
serial_results = benchmark_processing(test_selfies, use_parallel=False)
parallel_results = benchmark_processing(test_selfies, use_parallel=True)

import selfies as sf
from rdkit import Chem
from rdkit.Chem import Descriptors, Crippen

class SmartMolecularFilter:
    """Filter molecules based on drug-like properties."""
    
    def __init__(self):
        self.lipinski_rules = {
            'mw_max': 500,
            'logp_max': 5,
            'hbd_max': 5,
            'hba_max': 10
        }
    
    def check_lipinski(self, mol):
        """Check Lipinski's Rule of Five."""
        if not mol:
            return False
        
        mw = Descriptors.MolWt(mol)
        logp = Crippen.MolLogP(mol)
        hbd = Descriptors.NumHDonors(mol)
        hba = Descriptors.NumHAcceptors(mol)
        
        return (mw <= self.lipinski_rules['mw_max'] and
                logp <= self.lipinski_rules['logp_max'] and
                hbd <= self.lipinski_rules['hbd_max'] and
                hba <= self.lipinski_rules['hba_max'])
    
    def analyze_selfies(self, selfies_str):
        """Analyze a SELFIES string for drug-likeness."""
        try:
            smiles = sf.decoder(selfies_str)
            mol = Chem.MolFromSmiles(smiles)
            
            if not mol:
                return {'valid': False, 'reason': 'Invalid molecule'}
            
            # Calculate properties
            props = {
                'mw': Descriptors.MolWt(mol),
                'logp': Crippen.MolLogP(mol),
                'hbd': Descriptors.NumHDonors(mol),
                'hba': Descriptors.NumHAcceptors(mol),
                'tpsa': Descriptors.TPSA(mol),
                'rotatable_bonds': Descriptors.NumRotatableBonds(mol)
            }
            
            # Check filters
            lipinski_pass = self.check_lipinski(mol)
            
            return {
                'valid': True,
                'selfies': selfies_str,
                'smiles': smiles,
                'properties': props,
                'lipinski_compliant': lipinski_pass,
                'drug_like_score': self._calculate_drug_score(props)
            }
            
        except Exception as e:
            return {'valid': False, 'reason': f'Error: {str(e)}'}
    
    def _calculate_drug_score(self, props):
        """Simple drug-likeness score (0-1)."""
        score = 0
        if props['mw'] <= 500: score += 0.2
        if 0 <= props['logp'] <= 5: score += 0.2
        if props['hbd'] <= 5: score += 0.2
        if props['hba'] <= 10: score += 0.2
        if props['tpsa'] <= 140: score += 0.2
        return score

# Test the filter
filter_tool = SmartMolecularFilter()

# Test with some SELFIES
test_molecules = [
    "[C][C][Branch1][C][=O][O]",  # Acetic acid
    "[C][C][Branch1][C][=O][O][C][=C][C][=C][C][=C][Ring1][=Branch1][C][Branch1][C][=O][O]",  # Aspirin
    "[C][N][C][=O]"  # N-methylformamide
]

for selfies_str in test_molecules:
    result = filter_tool.analyze_selfies(selfies_str)
    if result['valid']:
        props = result['properties']
        print(f"SELFIES: {selfies_str}")
        print(f"SMILES:  {result['smiles']}")
        print(f"MW: {props['mw']:.1f}, LogP: {props['logp']:.2f}")
        print(f"Drug score: {result['drug_like_score']:.1f}/1.0")
        print(f"Lipinski: {'✓' if result['lipinski_compliant'] else '✗'}")
        print()

# Install dependencies
pip install selfies rdkit pillow

# Generate some random molecules
python -c "
import selfies as sf
import random
alphabet = list(sf.get_semantic_robust_alphabet())
for i in range(3):
    random_selfies = ''.join(random.choices(alphabet, k=8))
    smiles = sf.decoder(random_selfies)
    print(f'{i+1}. {random_selfies} → {smiles}')
"

# Visualize your favorite molecule (using the script below)
python selfies2png.py "[C][=C][C][=C][C][=C][Ring1][=Branch1]" -o benzene.png

#!/usr/bin/env python
# -*- coding: utf-8 -*-

"""
SELFIES to PNG Converter
=======================

A command-line utility to render SELFIES strings as 2D molecular images with
molecular formulas and proper subscript formatting.

This script demonstrates how to:
- Parse SELFIES strings and decode them to SMILES
- Generate 2D molecular coordinates using RDKit
- Create publication-quality molecular images
- Add custom legends with molecular formulas
- Handle font rendering and subscripts

Example Usage:
    python selfies2png.py "[C][C][O]"  # Ethanol
    python selfies2png.py "[C][C][Branch1][C][C][C][C]" -o isobutane.png  # Isobutane
    python selfies2png.py "[C][=C][C][=C][C][=C][Ring1][=Branch1]" --size 800  # Benzene, larger image

Author: Hunter Heidenreich
Website: https://hunterheidenreich.com
"""

import argparse
import hashlib
import sys
from pathlib import Path

# SELFIES import
import selfies as sf

# RDKit imports
from rdkit import Chem
from rdkit.Chem import Draw, rdDepictor, rdMolDescriptors

# PIL imports for image manipulation
from PIL import Image, ImageDraw, ImageFont

# Constants for image configuration
DEFAULT_IMAGE_SIZE = 500
LEGEND_HEIGHT_RATIO = 0.08  # Legend height as ratio of image size
LEGEND_Y_OFFSET_RATIO = 0.02  # Y offset as ratio of image size
LEGEND_X_OFFSET_RATIO = 0.02  # X offset as ratio of image size
SUBSCRIPT_Y_OFFSET_RATIO = 0.006  # Subscript offset as ratio of image size

# Font size ratios based on image size
REGULAR_FONT_RATIO = 0.028  # Regular font size as ratio of image size
SUBSCRIPT_FONT_RATIO = 0.02  # Subscript font size as ratio of image size

# Font paths for different operating systems
FONT_PATHS = [
    "/System/Library/Fonts/Arial.ttf",  # macOS
    "/usr/share/fonts/truetype/arial.ttf",  # Linux
    "C:/Windows/Fonts/arial.ttf",  # Windows
]


def _calculate_dynamic_sizes(image_size: int):
    """Calculate dynamic sizing values based on image size."""
    return {
        'legend_height': int(image_size * LEGEND_HEIGHT_RATIO),
        'legend_y_offset': int(image_size * LEGEND_Y_OFFSET_RATIO),
        'legend_x_offset': int(image_size * LEGEND_X_OFFSET_RATIO),
        'subscript_y_offset': int(image_size * SUBSCRIPT_Y_OFFSET_RATIO),
        'regular_font_size': int(image_size * REGULAR_FONT_RATIO),
        'subscript_font_size': int(image_size * SUBSCRIPT_FONT_RATIO),
    }


def _load_fonts(regular_size: int, subscript_size: int):
    """Load system fonts for text rendering, with fallback to default font."""
    font_regular = None
    font_small = None

    for font_path in FONT_PATHS:
        try:
            font_regular = ImageFont.truetype(font_path, regular_size)
            font_small = ImageFont.truetype(font_path, subscript_size)
            break
        except (OSError, IOError):
            continue

    if font_regular is None:
        font_regular = ImageFont.load_default()
        font_small = ImageFont.load_default()

    return font_regular, font_small


def create_molecule_image(
    mol: Chem.Mol, selfies_string: str, size: int = DEFAULT_IMAGE_SIZE
) -> Image.Image:
    """
    Creates a molecule image with a legend showing molecular formula and SELFIES string.

    Args:
        mol: RDKit molecule object (already validated)
        selfies_string: Original SELFIES string for legend display
        size: Image size in pixels (square image)

    Returns:
        PIL Image object with molecule structure and formatted legend
    """
    # Calculate dynamic sizes based on image size
    sizes = _calculate_dynamic_sizes(size)
    
    rdDepictor.Compute2DCoords(mol)
    molecular_formula = rdMolDescriptors.CalcMolFormula(mol)

    mol_img = Draw.MolToImage(mol, size=(size, size))
    if mol_img.mode != "RGBA":
        mol_img = mol_img.convert("RGBA")

    # Calculate required legend height based on content - be more generous
    legend_height = _calculate_legend_height(molecular_formula, selfies_string, sizes, size)
    
    # Add extra safety margin to ensure text doesn't get cut off
    legend_height = max(legend_height, int(size * 0.15))  # At least 15% of image height
    
    # Calculate if we need extra width for the SELFIES text
    temp_img = Image.new("RGBA", (1, 1), "white")
    temp_draw = ImageDraw.Draw(temp_img)
    font_regular, _ = _load_fonts(sizes['regular_font_size'], sizes['subscript_font_size'])
    
    selfies_label = "SELFIES: "
    label_width = int(temp_draw.textlength(selfies_label, font=font_regular))
    selfies_width = int(temp_draw.textlength(selfies_string, font=font_regular))
    total_selfies_width = label_width + selfies_width + (sizes['legend_x_offset'] * 2)
    
    # Make image wider if needed to accommodate the full SELFIES string
    final_width = max(size, total_selfies_width)
    
    total_height = size + legend_height
    final_img = Image.new("RGBA", (final_width, total_height), "white")
    
    # Center the molecule image if we made the canvas wider
    mol_x_offset = (final_width - size) // 2
    final_img.paste(mol_img, (mol_x_offset, 0))

    draw = ImageDraw.Draw(final_img)
    font_regular, font_small = _load_fonts(sizes['regular_font_size'], sizes['subscript_font_size'])

    _draw_molecular_formula(draw, molecular_formula, font_regular, font_small, sizes, size, mol_x_offset)
    _draw_selfies_legend(draw, selfies_string, font_regular, sizes, molecular_formula, size, final_width)

    return final_img


def _draw_molecular_formula(
    draw: ImageDraw.Draw, formula: str, font_regular, font_small, sizes: dict, image_size: int, mol_x_offset: int = 0
) -> int:
    """Draw molecular formula with proper subscript formatting."""
    y_pos = image_size + sizes['legend_y_offset']
    x_pos = sizes['legend_x_offset'] + mol_x_offset

    draw.text((x_pos, y_pos), "Formula: ", fill="black", font=font_regular)
    x_pos += int(draw.textlength("Formula: ", font=font_regular))

    for char in formula:
        if char.isdigit():
            draw.text(
                (x_pos, y_pos + sizes['subscript_y_offset']), char, fill="black", font=font_small
            )
            x_pos += int(draw.textlength(char, font=font_small))
        else:
            draw.text((x_pos, y_pos), char, fill="black", font=font_regular)
            x_pos += int(draw.textlength(char, font=font_regular))

    return x_pos


def _calculate_legend_height(formula: str, selfies: str, sizes: dict, image_size: int) -> int:
    """Calculate the required height for the legend based on content."""
    # Create a temporary draw object to measure text
    temp_img = Image.new("RGBA", (1, 1), "white")
    temp_draw = ImageDraw.Draw(temp_img)
    font_regular, _ = _load_fonts(sizes['regular_font_size'], sizes['subscript_font_size'])
    
    # Calculate if SELFIES will fit on the same line as formula
    formula_text = f"Formula: {formula}"
    separator = " | SELFIES: "
    total_prefix_width = int(temp_draw.textlength(formula_text, font=font_regular)) + int(temp_draw.textlength(separator, font=font_regular))
    available_width_same_line = image_size - sizes['legend_x_offset'] - total_prefix_width - sizes['legend_x_offset']
    selfies_width = int(temp_draw.textlength(selfies, font=font_regular))
    
    base_height = sizes['legend_y_offset'] * 2 + sizes['regular_font_size']
    line_height = sizes['regular_font_size'] + 6  # Line spacing to match drawing function
    
    if selfies_width <= available_width_same_line:
        # Single line layout - add extra padding
        return base_height + 15
    else:
        # SELFIES goes on new line - calculate width available for SELFIES line
        selfies_label_width = int(temp_draw.textlength("SELFIES: ", font=font_regular))
        available_width_new_line = image_size - sizes['legend_x_offset'] - selfies_label_width - sizes['legend_x_offset']
        
        if selfies_width <= available_width_new_line:
            # SELFIES fits on one new line
            return base_height + line_height + 15
        else:
            # SELFIES needs wrapping - parse into tokens for accurate calculation
            tokens = []
            current_token = ""
            
            for char in selfies:
                if char == '[':
                    if current_token:
                        tokens.append(current_token)
                        current_token = ""
                    current_token = char
                elif char == ']':
                    current_token += char
                    tokens.append(current_token)
                    current_token = ""
                else:
                    current_token += char
            
            if current_token:
                tokens.append(current_token)
            
            # Simulate line wrapping
            lines_needed = 1
            current_line_width = 0
            
            for token in tokens:
                token_width = int(temp_draw.textlength(token, font=font_regular))
                if current_line_width + token_width > available_width_new_line and current_line_width > 0:
                    lines_needed += 1
                    current_line_width = token_width
                else:
                    current_line_width += token_width
            
            # Add extra height for the second line (SELFIES line) plus wrapped lines plus extra padding
            total_height = base_height + line_height + (line_height * (lines_needed - 1)) + 25  # Extra generous padding
            return total_height


def _draw_wrapped_text(
    draw: ImageDraw.Draw, text: str, start_x: int, start_y: int, max_width: int, font, sizes: dict
) -> None:
    """Draw text with word wrapping at bracket boundaries for SELFIES strings."""
    # For SELFIES, we want to break at bracket boundaries to maintain readability
    # Split SELFIES into tokens (brackets and their contents)
    tokens = []
    current_token = ""
    
    for char in text:
        if char == '[':
            if current_token:
                tokens.append(current_token)
                current_token = ""
            current_token = char
        elif char == ']':
            current_token += char
            tokens.append(current_token)
            current_token = ""
        else:
            current_token += char
    
    if current_token:
        tokens.append(current_token)
    
    # Now draw tokens, wrapping as needed
    x_pos = start_x
    y_pos = start_y
    line_height = sizes['regular_font_size'] + 4  # More generous line spacing
    
    for token in tokens:
        token_width = int(draw.textlength(token, font=font))
        
        # If this token would exceed the line, wrap to next line
        if x_pos + token_width > start_x + max_width and x_pos > start_x:
            y_pos += line_height
            x_pos = start_x
        
        draw.text((x_pos, y_pos), token, fill="black", font=font)
        x_pos += token_width


def _draw_selfies_legend(
    draw: ImageDraw.Draw, selfies: str, font_regular, sizes: dict, formula: str, 
    original_mol_size: int, final_image_width: int
) -> None:
    """Add SELFIES string to the image legend with proper text wrapping."""
    y_pos = original_mol_size + sizes['legend_y_offset']
    
    # ALWAYS put SELFIES on its own line - no more trying to fit on same line as formula
    selfies_y_pos = y_pos + sizes['regular_font_size'] + 6
    selfies_label = "SELFIES: "
    label_width = int(draw.textlength(selfies_label, font=font_regular))
    
    draw.text((sizes['legend_x_offset'], selfies_y_pos), selfies_label, fill="black", font=font_regular)
    
    # Use wrapped text drawing for ALL SELFIES strings
    selfies_x_pos = sizes['legend_x_offset'] + label_width
    available_width_for_selfies = final_image_width - selfies_x_pos - sizes['legend_x_offset']
    _draw_wrapped_text(
        draw, selfies, selfies_x_pos, selfies_y_pos, 
        available_width_for_selfies, font_regular, sizes
    )


def selfies_to_png(
    selfies_string: str, output_file: str, size: int = DEFAULT_IMAGE_SIZE
) -> None:
    """
    Convert a SELFIES string to a PNG image with molecular formula legend.

    Args:
        selfies_string: Valid SELFIES string representing a molecule
        output_file: Path where the PNG image will be saved
        size: Square image dimensions in pixels

    Raises:
        ValueError: If SELFIES string is invalid or size is non-positive
        IOError: If file cannot be saved to the specified location
    """
    if not selfies_string or not selfies_string.strip():
        raise ValueError("SELFIES string cannot be empty")

    if size <= 0:
        raise ValueError(f"Image size must be positive, got: {size}")

    output_path = Path(output_file)
    output_path.parent.mkdir(parents=True, exist_ok=True)

    # Decode SELFIES to SMILES
    try:
        smiles_string = sf.decoder(selfies_string.strip())
    except Exception as e:
        raise ValueError(
            f"Invalid SELFIES string: '{selfies_string}'. "
            f"SELFIES decoding error: {e}"
        )

    if not smiles_string:
        raise ValueError(
            f"SELFIES string '{selfies_string}' decoded to empty SMILES. "
            f"Please check the SELFIES syntax."
        )

    mol = Chem.MolFromSmiles(smiles_string)
    if mol is None:
        raise ValueError(
            f"Invalid SELFIES string: '{selfies_string}'. "
            f"Decoded to SMILES '{smiles_string}' which is not valid. "
            f"Please check the syntax and try again."
        )

    img = create_molecule_image(mol, selfies_string.strip(), size)

    try:
        img.save(output_file, "PNG", optimize=True)
        print(f"Image successfully saved to: {output_file}")
    except Exception as e:
        raise IOError(f"Failed to save image to '{output_file}': {e}")


def create_safe_filename(selfies_string: str) -> str:
    """
    Generate a filesystem-safe filename from a SELFIES string using MD5 hash.

    Args:
        selfies_string: The input SELFIES string

    Returns:
        A safe filename ending with .png
    """
    clean_selfies = selfies_string.strip()
    hasher = hashlib.md5(clean_selfies.encode("utf-8"))
    return f"{hasher.hexdigest()}.png"


def main() -> None:
    """Command-line interface for the SELFIES to PNG converter."""
    parser = argparse.ArgumentParser(
        description="Convert SELFIES strings to publication-quality PNG images with molecular formulas.",
        epilog="""
Examples:
  %(prog)s "[C][C][O]"                                    # Ethanol with auto-generated filename
  %(prog)s "[C][C][Branch1][C][C][C][C]"                  # Isobutane with auto-generated filename
  %(prog)s "[C][=C][C][=C][C][=C][Ring1][=Branch1]" -o benzene.png   # Benzene with custom filename  
  %(prog)s "[C][C][O]" --size 800                         # Ethanol with larger image size

Common SELFIES patterns:
  [C][C][O]                                  - Ethanol
  [C][C][Branch1][C][=O][O]                  - Acetic acid
  [C][=C][C][=C][C][=C][Ring1][=Branch1]     - Benzene
  [C][C][Branch1][C][C][C]                   - Isobutane
  [N][C][Branch1][C][=O][C][=C][C][=C][C][=C][Ring1][=Branch1]  - Benzamide
        """,
        formatter_class=argparse.RawDescriptionHelpFormatter,
    )

    parser.add_argument(
        "selfies",
        type=str,
        help="SELFIES string of the molecule to visualize (e.g., '[C][C][O]' for ethanol)",
    )

    parser.add_argument(
        "-o",
        "--output",
        type=str,
        metavar="FILE",
        help="Output PNG filename. If not provided, generates a unique filename "
        "based on the SELFIES string hash. Extension .png will be added if missing.",
    )

    parser.add_argument(
        "-s",
        "--size",
        type=int,
        default=DEFAULT_IMAGE_SIZE,
        metavar="PIXELS",
        help=f"Square image size in pixels (default: {DEFAULT_IMAGE_SIZE}). "
        f"Typical values: 300 (small), 500 (medium), 800 (large).",
    )

    parser.add_argument(
        "-v",
        "--verbose",
        action="store_true",
        help="Enable verbose output for debugging",
    )

    args = parser.parse_args()

    if args.verbose:
        print(f"Input SELFIES: {args.selfies}")
        print(f"Image size: {args.size}x{args.size} pixels")

    if args.output:
        output_filename = (
            args.output
            if args.output.lower().endswith(".png")
            else f"{args.output}.png"
        )
        if args.verbose:
            print(f"Using custom filename: {output_filename}")
    else:
        output_filename = create_safe_filename(args.selfies)
        if args.verbose:
            print(f"Generated filename: {output_filename}")

    try:
        selfies_to_png(args.selfies, output_filename, args.size)

        if args.verbose:
            # Decode and show the SMILES for reference
            try:
                decoded_smiles = sf.decoder(args.selfies.strip())
                print(f"Decoded SMILES: {decoded_smiles}")
            except Exception:
                pass
            print("Conversion completed successfully!")

    except ValueError as e:
        print(f"Input Error: {e}", file=sys.stderr)
        print("Tip: Check your SELFIES string syntax", file=sys.stderr)
        sys.exit(1)

    except IOError as e:
        print(f"File Error: {e}", file=sys.stderr)
        print("Tip: Check file permissions and disk space", file=sys.stderr)
        sys.exit(2)

    except ImportError as e:
        print(f"Dependencies Error: {e}", file=sys.stderr)
        print(
            "Tip: Install required packages with 'pip install rdkit pillow selfies'",
            file=sys.stderr,
        )
        sys.exit(3)

    except Exception as e:
        print(f"Unexpected Error: {e}", file=sys.stderr)
        print("Tip: Please report this issue if it persists", file=sys.stderr)
        sys.exit(4)


if __name__ == "__main__":
    main()

# Basic usage - generates auto-named file
python selfies2png.py "[C][C][O]"

# Custom filename and size
python selfies2png.py "[C][=C][C][=C][C][=C][Ring1][=Branch1]" -o benzene.png --size 800

# Complex molecules with branching
python selfies2png.py "[C][C][Branch1][C][=O][O][C][=C][C][=C][C][=C][Ring1][=Branch1][C][Branch1][C][=O][O]" -o aspirin.png