@@ -3539,30 +3539,23 @@ tsk_treeseq_site_allele_frequency_spectrum(const tsk_treeseq_t *self,
35393539 return ret ;
35403540}
35413541
3542- static int TSK_WARN_UNUSED
3542+ static void
35433543tsk_treeseq_update_branch_afs (const tsk_treeseq_t * self , tsk_id_t u , double right ,
35443544 double * restrict last_update , const double * restrict time , tsk_id_t * restrict parent ,
3545- const double * counts , tsk_size_t num_sample_sets , tsk_size_t num_time_windows ,
3546- const double * time_windows , tsk_size_t window_index , const tsk_size_t * result_dims ,
3547- tsk_flags_t options , double * result )
3545+ tsk_size_t * restrict coordinate , const double * counts , tsk_size_t num_sample_sets ,
3546+ tsk_size_t num_time_windows , const double * time_windows , tsk_size_t window_index ,
3547+ const tsk_size_t * result_dims , tsk_flags_t options , double * result )
35483548{
3549- int ret = 0 ;
35503549 tsk_size_t afs_size ;
35513550 tsk_size_t k ;
35523551 tsk_size_t time_window_index ;
35533552 double * afs ;
3554- // note: moving this malloc outside this function doesn't speed things up
3555- tsk_size_t * coordinate = tsk_malloc (num_sample_sets * sizeof (* coordinate ));
35563553 bool polarised = !!(options & TSK_STAT_POLARISED );
35573554 const double * count_row = GET_2D_ROW (counts , num_sample_sets + 1 , u );
35583555 double x = 0 ;
35593556 double t_u , t_v ;
35603557 double tw_branch_length = 0 ;
35613558 const tsk_size_t all_samples = (tsk_size_t ) count_row [num_sample_sets ];
3562- if (coordinate == NULL ) {
3563- ret = tsk_trace_error (TSK_ERR_NO_MEMORY );
3564- goto out ;
3565- }
35663559 if (parent [u ] != TSK_NULL ) {
35673560 t_u = time [u ];
35683561 t_v = time [parent [u ]];
@@ -3589,9 +3582,6 @@ tsk_treeseq_update_branch_afs(const tsk_treeseq_t *self, tsk_id_t u, double righ
35893582 }
35903583 }
35913584 last_update [u ] = right ;
3592- out :
3593- tsk_safe_free (coordinate );
3594- return ret ;
35953585}
35963586
35973587static int
@@ -3616,6 +3606,7 @@ tsk_treeseq_branch_allele_frequency_spectrum(const tsk_treeseq_t *self,
36163606 tsk_id_t * restrict parent = tsk_malloc (num_nodes * sizeof (* parent ));
36173607 double * restrict last_update = tsk_calloc (num_nodes , sizeof (* last_update ));
36183608 double * restrict branch_length = tsk_calloc (num_nodes , sizeof (* branch_length ));
3609+ tsk_size_t * restrict coordinate = tsk_malloc (num_sample_sets * sizeof (* coordinate ));
36193610 tsk_id_t tj , tk , h ;
36203611 double t_left , t_right , w_right ;
36213612 const tsk_size_t K = num_sample_sets + 1 ;
@@ -3625,7 +3616,7 @@ tsk_treeseq_branch_allele_frequency_spectrum(const tsk_treeseq_t *self,
36253616 goto out ;
36263617 }
36273618
3628- if (parent == NULL || last_update == NULL ) {
3619+ if (parent == NULL || last_update == NULL || coordinate == NULL ) {
36293620 ret = tsk_trace_error (TSK_ERR_NO_MEMORY );
36303621 goto out ;
36313622 }
@@ -3643,19 +3634,13 @@ tsk_treeseq_branch_allele_frequency_spectrum(const tsk_treeseq_t *self,
36433634 tk ++ ;
36443635 u = edge_child [h ];
36453636 v = edge_parent [h ];
3646- ret = tsk_treeseq_update_branch_afs (self , u , t_left , last_update , node_time ,
3647- parent , counts , num_sample_sets , num_time_windows , time_windows ,
3648- window_index , result_dims , options , result );
3649- if (ret != 0 ) {
3650- goto out ;
3651- }
3637+ tsk_treeseq_update_branch_afs (self , u , t_left , last_update , node_time ,
3638+ parent , coordinate , counts , num_sample_sets , num_time_windows ,
3639+ time_windows , window_index , result_dims , options , result );
36523640 while (v != TSK_NULL ) {
3653- ret = tsk_treeseq_update_branch_afs (self , v , t_left , last_update ,
3654- node_time , parent , counts , num_sample_sets , num_time_windows ,
3641+ tsk_treeseq_update_branch_afs (self , v , t_left , last_update , node_time ,
3642+ parent , coordinate , counts , num_sample_sets , num_time_windows ,
36553643 time_windows , window_index , result_dims , options , result );
3656- if (ret != 0 ) {
3657- goto out ;
3658- }
36593644 update_state (counts , K , v , u , -1 );
36603645 v = parent [v ];
36613646 }
@@ -3671,12 +3656,9 @@ tsk_treeseq_branch_allele_frequency_spectrum(const tsk_treeseq_t *self,
36713656 parent [u ] = v ;
36723657 branch_length [u ] = node_time [v ] - node_time [u ];
36733658 while (v != TSK_NULL ) {
3674- ret = tsk_treeseq_update_branch_afs (self , v , t_left , last_update ,
3675- node_time , parent , counts , num_sample_sets , num_time_windows ,
3659+ tsk_treeseq_update_branch_afs (self , v , t_left , last_update , node_time ,
3660+ parent , coordinate , counts , num_sample_sets , num_time_windows ,
36763661 time_windows , window_index , result_dims , options , result );
3677- if (ret != 0 ) {
3678- goto out ;
3679- }
36803662 update_state (counts , K , v , u , +1 );
36813663 v = parent [v ];
36823664 }
@@ -3695,12 +3677,9 @@ tsk_treeseq_branch_allele_frequency_spectrum(const tsk_treeseq_t *self,
36953677 /* Flush the contributions of all nodes to the current window */
36963678 for (u = 0 ; u < (tsk_id_t ) num_nodes ; u ++ ) {
36973679 tsk_bug_assert (last_update [u ] < w_right );
3698- ret = tsk_treeseq_update_branch_afs (self , u , w_right , last_update ,
3699- node_time , parent , counts , num_sample_sets , num_time_windows ,
3680+ tsk_treeseq_update_branch_afs (self , u , w_right , last_update , node_time ,
3681+ parent , coordinate , counts , num_sample_sets , num_time_windows ,
37003682 time_windows , window_index , result_dims , options , result );
3701- if (ret != 0 ) {
3702- goto out ;
3703- }
37043683 }
37053684 window_index ++ ;
37063685 }
0 commit comments