Skip to content

Commit

Permalink
SIRENA update (v.8.0.1): Solved bug (freeing 'model' in runEnergy)
Browse files Browse the repository at this point in the history
  • Loading branch information
bcobo committed Nov 6, 2023
1 parent 41c7311 commit 938813b
Show file tree
Hide file tree
Showing 5 changed files with 87 additions and 18 deletions.
74 changes: 71 additions & 3 deletions libsixt/initSIRENA.c
Expand Up @@ -54,6 +54,7 @@ int checkXmls(struct Parameters* const par)
}
else // or read full Primary HDU and store it in 'libheaderPrimary' and calculate XML checksum
{
status = EXIT_SUCCESS;
int numberkeywords;
char *libheaderPrimary = NULL;
if (fits_hdr2str(libptr, 0, NULL, 0,&libheaderPrimary, &numberkeywords, &status))
Expand Down Expand Up @@ -453,7 +454,7 @@ int getSamplingrate_trigreclength (char* inputFile, struct Parameters par, doubl


/***** SECTION 6 ************************************************************
* fillReconstructInitSIRENA: This function reads the grading data from the XML file
* fillReconstructInitSIRENAGrading: This function reads the grading data from the XML file
* and store it in 'reconstruct_init_sirena->grading'
*
* Parameters:
Expand Down Expand Up @@ -538,7 +539,7 @@ int fillReconstructInitSIRENAGrading (struct Parameters par, AdvDet *det, Recons
}
}

int OFlengthvsposti = 0;
/*int OFlengthvsposti = 0;
if ((par.preBuffer == 1) && (par.opmode == 1))
{
for (int i=0;i<(*reconstruct_init_sirena)->grading->ngrades;i++)
Expand All @@ -550,11 +551,14 @@ int fillReconstructInitSIRENAGrading (struct Parameters par, AdvDet *det, Recons
}
}
if (OFlengthvsposti == 0)
printf("%s %d %s","par.: ",par.Pul,"\n");
printf("%s %d %s","par.OFLength: ",par.OFLength,"\n");
if ((OFlengthvsposti == 0) && ((*reconstruct_init_sirena)->pulse_length >= (*reconstruct_init_sirena)->OFLength)) // Not 0-padding
{
SIXT_ERROR("The grading/preBuffer info of the XML file does not match the OFLength input parameter");
return(EXIT_FAILURE);
}
}
}*/

return(status);
}
Expand Down Expand Up @@ -618,13 +622,77 @@ int callSIRENA_Filei(char* inputFile, SixtStdKeywords* keywords, ReconstructInit
if (status != EXIT_SUCCESS) return(EXIT_FAILURE);

// Iterate of records and run SIRENA
/*fitsfile *fptr;
//if (fits_open_file(&fptr, inputFile, READONLY, &status)) {
// fits_report_error(stderr, status);
// return status;
//}
fits_open_file(&fptr, inputFile, READONLY, &status);
fits_movnam_hdu(fptr, ANY_HDU,"TESRECORDS", 0, &status);
//if (status != EXIT_SUCCESS) return(EXIT_FAILURE);
int numrecords;
fits_read_key(fptr, TINT, "NAXIS2", &numrecords, NULL, &status);
fits_close_file(fptr, &status);
//printf("%s","Paso5\n");
//printf("%s %d %s","NAXIS2=",numrecords,"\n");
//printf("%s","Paso6\n");
//int i;
//int total=numrecords;
// for (i = 0; i <= total; i++) {
// float progress = (float)i / total;
// int bar_width = 50;
// printf("Simulating |");
// int pos = bar_width * progress;
// for (int j = 0; j < bar_width; j++) {
// if (j < pos)
// printf("=");
// else if (j == pos)
// printf(">");
// else
// printf(" ");
// }
// printf("| %.2f%%\r", progress * 100);
// fflush(stdout);
// usleep(100000); // Sleep for a short time to simulate progress
//}
//printf("\n");*/

fits_movnam_hdu(outfile->fptr, ANY_HDU,"EVENTS", 0, &status);
if (status != EXIT_SUCCESS) return(EXIT_FAILURE);

int lastRecord = 0, nrecord = 0, startRecordGroup = 0; //last record required for SIRENA library creation
// Build up TesEventList to recover the results of SIRENA
TesEventList* event_list = newTesEventListSIRENA(&status);
allocateTesEventListTrigger(event_list,par.EventListSize,&status);
if (status != EXIT_SUCCESS) return(EXIT_FAILURE);
while(getNextRecord(record_file,record,&lastRecord,&startRecordGroup,&status))
{
/*float progress = (float)nrecord / numrecords;
int bar_width = 50;
printf("Simulating |");
int pos = bar_width * progress;
for (int j = 0; j < bar_width; j++) {
if (j < pos)
printf("=");
else if (j == pos)
printf(">");
else
printf(" ");
}
if (nrecord == numrecords-1)
{
printf("| %.2f%%\r", progress * 100);
fflush(stdout);
usleep(100000); // Sleep for a short time to simulate progress
printf("\n");
}*/

nrecord = nrecord + 1;
/*if(nrecord < 5887)
{
Expand Down
4 changes: 2 additions & 2 deletions libsixt/integraSIRENA.cpp
Expand Up @@ -200,7 +200,7 @@
}
largeFilter = pulse_length;
}

reconstruct_init->library_collection = getLibraryCollection(library_file, opmode, hduPRECALWN, hduPRCLOFWM, largeFilter, filter_domain, pulse_length, energy_method, ofnoise, filter_method, oflib, &ofinterp, filtEev, lagsornot, reconstruct_init->preBuffer, pBi, posti, status);
reconstruct_init->library_collection->margin = 0.25; // (%) Margin to be applied when several energies in the library to choose the proper filter
if (*status)
Expand Down Expand Up @@ -2748,7 +2748,7 @@
if (ntemplates == 1) nOFs = nOFs-1; // -1 because the ENERGYcolumn
else nOFs = (nOFs-1)/2; // /2 because the AB column
}

if (nOFs == 0)
{
EP_PRINT_ERROR("The library has no fixed optimal filters",EPFAIL);
Expand Down
21 changes: 11 additions & 10 deletions libsixt/tasksSIRENA.cpp
Expand Up @@ -8364,6 +8364,7 @@ void runEnergy(TesRecord* record, int lastRecord, int nrecord, int trig_reclengt
int resize_mfvsposti = 0;

int extraSizeDueToLags = 0;
model = gsl_vector_alloc((*reconstruct_init)->pulse_length);
for (int i=0; i<(*pulsesInRecord)->ndetpulses ;i++)
{
log_debug("Pulse................................................ %d",i+1);
Expand All @@ -8381,7 +8382,6 @@ void runEnergy(TesRecord* record, int lastRecord, int nrecord, int trig_reclengt
message = "Cannot run routine pulseGrading";
EP_EXIT_ERROR(message,EPFAIL);
}
model = gsl_vector_alloc((*reconstruct_init)->pulse_length);

if (preBuffer == 1)
{
Expand Down Expand Up @@ -9167,8 +9167,8 @@ void runEnergy(TesRecord* record, int lastRecord, int nrecord, int trig_reclengt
} // End for
log_debug("After FOR");

gsl_vector_free(recordAux); recordAux = 0;
gsl_vector_free(model); model = 0;
if (recordAux != NULL) {gsl_vector_free(recordAux); recordAux = 0;}
if (model != NULL) {gsl_vector_free(model); model = 0;}

if (pulse_lowres != NULL) {gsl_vector_free(pulse_lowres); pulse_lowres = 0;}
if (filtergsl_lowres!= NULL) {gsl_vector_free(filtergsl_lowres); filtergsl_lowres = 0;}
Expand Down Expand Up @@ -9213,7 +9213,7 @@ void runEnergy(TesRecord* record, int lastRecord, int nrecord, int trig_reclengt


/***** SECTION BB ************************************************************
* th_runEenergy: Run energy calculation only in multithread mode
* th_runEnergy: Run energy calculation only in multithread mode
*****************************************************************************/
void th_runEnergy(TesRecord* record, int nrecord, int trig_reclength,
ReconstructInitSIRENA** reconstruct_init,
Expand Down Expand Up @@ -9415,6 +9415,7 @@ void th_runEnergy(TesRecord* record, int nrecord, int trig_reclength,
int resize_mfvsposti = 0;

int extraSizeDueToLags = 0;
model =gsl_vector_alloc((*reconstruct_init)->pulse_length);
for (int i=0; i<(*pulsesInRecord)->ndetpulses ;i++)
{
tstartSamplesRecord = (*pulsesInRecord)->pulses_detected[i].TstartSamples;
Expand All @@ -9431,7 +9432,6 @@ void th_runEnergy(TesRecord* record, int nrecord, int trig_reclength,
message = "Cannot run routine pulseGrading";
EP_EXIT_ERROR(message,EPFAIL);
}
model =gsl_vector_alloc((*reconstruct_init)->pulse_length);

if (preBuffer == 1)
{
Expand Down Expand Up @@ -10206,11 +10206,10 @@ void th_runEnergy(TesRecord* record, int nrecord, int trig_reclength,
(*pulsesInRecord)->pulses_detected[i].phi = -999.0;
(*pulsesInRecord)->pulses_detected[i].lagsShift = -999.0;
}

gsl_vector_free(model); model = 0;
} // End for

gsl_vector_free(recordAux); recordAux = 0;
if (recordAux != NULL) {gsl_vector_free(recordAux); recordAux = 0;}
if (model != NULL) {gsl_vector_free(model); model = 0;}

gsl_vector_free(pulse_lowres); pulse_lowres = 0;
gsl_vector_free(filtergsl_lowres); filtergsl_lowres = 0;
Expand Down Expand Up @@ -11822,6 +11821,7 @@ int pulseGrading (ReconstructInitSIRENA *reconstruct_init, int tstart, int grade
break;
}
}

if (reconstruct_init->preBuffer == 1)
{
if (tstart < gsl_matrix_get(reconstruct_init->grading->gradeData,*pulseGrade-1,2))
Expand Down Expand Up @@ -11954,8 +11954,8 @@ int pulseGrading (ReconstructInitSIRENA *reconstruct_init, int tstart, int grade
//if (((grade2 < gradelim_pre) || (grade1 == -1)) && (OFlength_strategy != 2)) *pulseGrade = -1;
if (grade2 < gradelim_pre) *pulseGrade = -1;

log_debug("*pulseGrade %d", *pulseGrade);
log_debug("*OFLength %d", *OFlength);
log_debug("pulseGrading *pulseGrade %d", *pulseGrade);
log_debug("pulseGrading *OFLength %d", *OFlength);

message.clear();

Expand Down Expand Up @@ -12093,6 +12093,7 @@ int calculateEnergy (gsl_vector *pulse, gsl_vector *filter, gsl_vector_complex *
log_debug("pulse->size: %i",pulse->size);
log_debug("productSize: %i",productSize);
}
//if (LowRes!=1){for (int i=0;i<10;i++) cout<<i<<" "<<gsl_vector_get(filter,i)<<endl;}

gsl_vector *vector;

Expand Down
4 changes: 2 additions & 2 deletions libsixt/versionSIRENA.h
Expand Up @@ -22,11 +22,11 @@
// CANTABRIA (CSIC-UC) with funding from the Spanish Ministry of Science and
// Innovation (MICINN)
//
// DATE: 2023/10/10, 09:52:56
// DATE: 2023/11/06, 12:19:50

#ifndef SIRENA_VERSION_H
#define SIRENA_VERSION_H

#define SIRENA_VERSION "8.0.0"
#define SIRENA_VERSION "8.0.1"

#endif
2 changes: 1 addition & 1 deletion tools/teslib/teslib.c
Expand Up @@ -91,7 +91,7 @@
int teslib_main() {
printf("Running TESLIB\n");
time_t ttstart = time(0);

// Containing all programm parameters read by PIL.
struct Parameters par;
par.hduPRCLOFWM = 0; // Debugger complains about an initialized variable (only the boolean type)
Expand Down

0 comments on commit 938813b

Please sign in to comment.