/* sim_ann.h contains the structures and definitions for simulated
   annealing optimisation */
/* RJ 12 September 2000 */

#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>
#include <malloc.h>
#include <time.h>
#include <memory.h>

#define SUCCESS 1
#define FAILURE 0

#define Ndvh_bins 1000  // number of bins to use for the calculation of all dvhs

#define SMARTBEAMLETS 1  // if "1", beamlet voxels that are not used in the optimization
                         // are not loaded, but the full beamlets must be reloaded every
						 // time the dose distribution is written to file.

typedef struct
{
  int x_count;
  int y_count;
  int z_count;
  int Nvox;
  float *matrix;
} GRID;

typedef struct
{
  int size;
  int *indices;
  float *data;
} SPARSE_ARRAY;

typedef struct
{ 
	int num_beamlets;
	int num_beamlet_batches;
    int x_count;
    int y_count;
    int z_count;
	int Nvox;  // number of voxels
    float *beam_weight;
    float *initial_beam_weight;
	SPARSE_ARRAY *array;
} GRID_ARRAY;

typedef struct
{
	int Nbins;         // number of dose bins
	float del_d;       // spacing between dose bins
	float dmax;        // maximum dose bin boundary in the dvh
	float del_d_plus;  // dose between the dVPlus threshold and the actual DVH
	float del_d_minus; // dose between the dVMinus threshold and the actual DVH
	float dvh[Ndvh_bins];        // dvh itself
} DVH;

typedef struct
{
	char name[100];

	float alpha;   // tissue importance factor

	DVH dvh;                    // dvh structure
	// DVH overdose penalty factors
	float betaVPlus;             // importance of dvh plus penalty term
	float dVPlus;
	float vPlus;
	// DVH underdose penalty factors
	float betaVMinus;            // importance of dvh minus penalty term
	float dVMinus;
	float vMinus;

	char dosePlusFilename[200]; // name of the file that contains overdose penalty information
	char doseMinusFilename[200]; // name of the file that contains underdose penalty information

	// mask indices
	int Nind;         // number of indices in the tissue mask
	int *ind;         // indices in which there is a non-zero dose prescription
	float betaPlus;   // overdose penalty weighting factor
	float *dPlus;     // overdose penalty dose threshold
	float betaMinus;  // underdose penalty weighting factor
	float *dMinus;    // underdose penalty weighting factor
} TISSUE;

typedef struct
{
	int tissueInd;  // index of prescription tissue being represented
	float dPlus;    // overdose penalty threshold
	float dMinus;   // underdose penalty threshold
} PRESCNODE;

typedef struct
{
	// grid dimensions
	int x_count;
	int y_count;
	int z_count;
	int Nvox;   // total number of voxels (just x_count*y_count*z_count)
	
	// the two following grids have the above dimensions
	PRESCNODE **array; // grid full of prescnode arrays
	int *numNodes;     // number of prescnodes at each point in grid
} PRESCGRID;

typedef struct
{
	int Niterations;
	int Nperbatch;
	float modFactor;
	char prescription_filename[200];
	char initial_beam_weights_filename[200];
	char beamlet_header_file[200];
	char dose_batch_base_name[200];
	char dose_batch_extension[200];
	char weight_batch_base_name[200];
	char weight_batch_extension[200];
	char obj_func_name[200];
	char bix_dir[200];  // directory containing beamlet batches
	char bix_batch_base_filename[200];
	char bix_batch_extension[200];
}  OPT_PARAM;   // optimization parameters

typedef struct
{
	int Ntissue;   // number of tissue types
	int x_count;   // dimensions of the tissue grids
	int y_count;
	int z_count;
	int Nvox;

	TISSUE *array; // array of tissue prescription structures
}  PRESC;

// function prototypes
int read_in_opt_parameters(OPT_PARAM *, char *);
int read_in_beamlets(GRID_ARRAY *, char *, char *, int *, OPT_PARAM *, GRID *);
int read_in_prescription(PRESC *, char *);
int calculate_dose_diff(GRID_ARRAY *, float, int, GRID *);
int calculate_dose(GRID_ARRAY *, GRID *);
int output_dose_weights(GRID_ARRAY *, GRID *, char *, char *);
int calc_update_factor(GRID_ARRAY *, GRID *, PRESC *, int, float *,PRESCGRID *);
int calc_dvh(GRID *, int, int *, DVH *);
float calc_obj_func(GRID_ARRAY *, GRID *, PRESC *);
int compare(void *, void *);
int calcAllDvhs(GRID *, PRESC *);
int find_used_voxels(PRESC *, int *);
int fillPrescriptionGrid(PRESC *, PRESCGRID *);
int loadBeamletsCalcDose(OPT_PARAM *, GRID *, GRID_ARRAY *);
int truncateBeamWeights(float *, int, float);