aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--src/likelihood.c9
-rw-r--r--src/likelihood.h7
-rw-r--r--src/process_hkl.c37
-rw-r--r--src/utils.c90
-rw-r--r--src/utils.h22
5 files changed, 144 insertions, 21 deletions
diff --git a/src/likelihood.c b/src/likelihood.c
index 95f005c7..98ffcba0 100644
--- a/src/likelihood.c
+++ b/src/likelihood.c
@@ -19,17 +19,20 @@
void detwin_intensities(const double *model, double *new_pattern,
const unsigned int *model_counts,
- unsigned int *new_counts)
+ ReflItemList *items)
{
/* Placeholder... */
}
void scale_intensities(const double *model, double *new_pattern,
const unsigned int *model_counts,
- unsigned int *new_counts, double f0, int f0_valid)
+ ReflItemList *items, double f0, int f0_valid)
{
double s;
unsigned int i;
+ unsigned int *new_counts;
+
+ new_counts = items_to_counts(items);
s = stat_scale_intensity(model, model_counts, new_pattern, new_counts);
if ( f0_valid ) printf("%f %f\n", s, f0);
@@ -41,4 +44,6 @@ void scale_intensities(const double *model, double *new_pattern,
for ( i=0; i<LIST_SIZE; i++ ) {
new_counts[i] *= s;
}
+
+ free(new_counts);
}
diff --git a/src/likelihood.h b/src/likelihood.h
index 0a2bb407..ad6585eb 100644
--- a/src/likelihood.h
+++ b/src/likelihood.h
@@ -17,13 +17,16 @@
#define LIKELIHOOD_H
+#include "utils.h"
+
+
extern void detwin_intensities(const double *model, double *new_pattern,
const unsigned int *model_counts,
- unsigned int *new_counts);
+ ReflItemList *items);
extern void scale_intensities(const double *model, double *new_pattern,
const unsigned int *model_counts,
- unsigned int *new_counts, double f0,
+ ReflItemList *items, double f0,
int f0_valid);
diff --git a/src/process_hkl.c b/src/process_hkl.c
index 77d5056a..deeb77f5 100644
--- a/src/process_hkl.c
+++ b/src/process_hkl.c
@@ -188,17 +188,21 @@ static void process_reflections(double *ref, unsigned int *counts,
static void merge_pattern(double *model, const double *new,
unsigned int *model_counts,
- const unsigned int *counts, int mo, int sum)
+ ReflItemList *items, int mo, int sum)
{
- signed int h, k, l;
+ int i;
- for ( l=-INDMAX; l<INDMAX; l++ ) {
- for ( k=-INDMAX; k<INDMAX; k++ ) {
- for ( h=-INDMAX; h<INDMAX; h++ ) {
+ for ( i=0; i<num_items(items); i++ ) {
double intensity;
+ signed int h, k, l;
+ struct refl_item *item;
- if ( lookup_count(counts, h, k, l) == 0 ) continue;
+ item = get_item(items, i);
+
+ h = item->h;
+ k = item->k;
+ l = item->l;
intensity = lookup_intensity(new, h, k, l);
@@ -217,8 +221,6 @@ static void merge_pattern(double *model, const double *new,
}
}
- }
- }
}
@@ -243,7 +245,6 @@ int main(int argc, char *argv[])
int config_scale = 0;
char *intfile = NULL;
double *new_pattern = NULL;
- unsigned int *new_counts = NULL;
unsigned int n_total_patterns;
unsigned int *truecounts = NULL;
char *sym = NULL;
@@ -251,6 +252,7 @@ int main(int argc, char *argv[])
float f0;
int f0_valid;
int n_nof0 = 0;
+ ReflItemList *items;
/* Long options */
const struct option longopts[] = {
@@ -344,7 +346,7 @@ int main(int argc, char *argv[])
model_counts = new_list_count();
cell = load_cell_from_pdb(pdb);
new_pattern = new_list_intensity();
- new_counts = new_list_count();
+ items = new_items();
if ( strcmp(filename, "-") == 0 ) {
fh = stdin;
@@ -393,7 +395,7 @@ int main(int argc, char *argv[])
/* Detwin before scaling */
if ( config_detwin ) {
detwin_intensities(model, new_pattern,
- model_counts, new_counts);
+ model_counts, items);
}
/* Assume a default I0 if we don't have one by now */
@@ -405,13 +407,13 @@ int main(int argc, char *argv[])
/* Scale if requested */
if ( config_scale ) {
scale_intensities(model, new_pattern,
- model_counts, new_counts, f0,
+ model_counts, items, f0,
f0_valid);
}
/* Start of second or later pattern */
merge_pattern(model, new_pattern, model_counts,
- new_counts, config_maxonly, config_sum);
+ items, config_maxonly, config_sum);
if ( (trueref != NULL) && config_every
&& (n_patterns % config_every == 0) ) {
@@ -424,9 +426,8 @@ int main(int argc, char *argv[])
if ( n_patterns == config_stopafter ) break;
- zero_list_count(new_counts);
-
n_patterns++;
+ clear_items(items);
progress_bar(n_patterns, n_total_patterns, "Merging");
@@ -449,12 +450,12 @@ int main(int argc, char *argv[])
if ( (h==0) && (k==0) && (l==0) ) continue;
- if ( lookup_count(new_counts, h, k, l) != 0 ) {
+ if ( find_item(items, h, k, l) != 0 ) {
ERROR("More than one measurement for %i %i %i in"
" pattern number %i\n", h, k, l, n_patterns);
}
set_intensity(new_pattern, h, k, l, intensity);
- set_count(new_counts, h, k, l, 1);
+ add_item(items, h, k, l);
} while ( rval != NULL );
@@ -480,5 +481,7 @@ int main(int argc, char *argv[])
STATUS("There were %u patterns.\n", n_patterns);
STATUS("%i had no f0 valid value.\n", n_nof0);
+ delete_items(items);
+
return 0;
}
diff --git a/src/utils.c b/src/utils.c
index e3b8b667..de820cf3 100644
--- a/src/utils.c
+++ b/src/utils.c
@@ -283,3 +283,93 @@ int assplode(const char *a, const char *delims, char ***pbits,
*pbits = bits;
return n;
}
+
+
+struct _reflitemlist {
+ struct refl_item *items;
+ int n_items;
+ int max_items;
+};
+
+void clear_items(ReflItemList *items)
+{
+ items->n_items = 0;
+}
+
+static void alloc_items(ReflItemList *items)
+{
+ items->items = realloc(items->items,
+ items->max_items*sizeof(struct refl_item));
+}
+
+ReflItemList *new_items()
+{
+ ReflItemList *new;
+ new = malloc(sizeof(ReflItemList));
+ new->max_items = 1024;
+ new->n_items = 0;
+ new->items = NULL;
+ alloc_items(new);
+ return new;
+}
+
+void delete_items(ReflItemList *items)
+{
+ if ( items->items != NULL ) free(items->items);
+ free(items);
+}
+
+void add_item(ReflItemList *items,
+ signed int h, signed int k, signed int l)
+{
+ if ( items->n_items == items->max_items ) {
+ items->max_items += 1024;
+ alloc_items(items);
+ }
+
+ items->items[items->n_items].h = h;
+ items->items[items->n_items].k = k;
+ items->items[items->n_items].l = l;
+ items->n_items++;
+}
+
+int find_item(ReflItemList *items,
+ signed int h, signed int k, signed int l)
+{
+ int i;
+
+ for ( i=0; i<items->n_items; i++ ) {
+ if ( items->items[i].h != h ) continue;
+ if ( items->items[i].k != k ) continue;
+ if ( items->items[i].l != l ) continue;
+ return 1;
+ }
+ return 0;
+}
+
+struct refl_item *get_item(ReflItemList *items, int i)
+{
+ if ( i >= items->n_items ) return NULL;
+ return &items->items[i];
+}
+
+int num_items(ReflItemList *items)
+{
+ return items->n_items;
+}
+
+unsigned int *items_to_counts(ReflItemList *items)
+{
+ unsigned int *c;
+ int i;
+
+ c = new_list_count();
+
+ for ( i=0; i<num_items(items); i++ ) {
+ struct refl_item *r;
+ r = get_item(items, i);
+ set_count(c, r->h, r->k, r->l, 1);
+ }
+
+ return c;
+}
diff --git a/src/utils.h b/src/utils.h
index b01725dd..14ef45a7 100644
--- a/src/utils.h
+++ b/src/utils.h
@@ -176,6 +176,28 @@ static inline double angle_between(double x1, double y1, double z1,
#include "list_tmp.h"
+/* ----------- Reflection lists indexed by sequence (not indices) ----------- */
+
+typedef struct _reflitemlist ReflItemList; /* Opaque */
+
+struct refl_item {
+ signed int h;
+ signed int k;
+ signed int l;
+};
+
+extern void clear_items(ReflItemList *items);
+extern ReflItemList *new_items(void);
+extern void delete_items(ReflItemList *items);
+extern void add_item(ReflItemList *items,
+ signed int h, signed int k, signed int l);
+extern int find_item(ReflItemList *items,
+ signed int h, signed int k, signed int l);
+extern struct refl_item *get_item(ReflItemList *items, int i);
+extern int num_items(ReflItemList *items);
+extern unsigned int *items_to_counts(ReflItemList *items);
+
+
/* ------------------------------ Message macros ---------------------------- */
#define ERROR(...) fprintf(stderr, __VA_ARGS__)