p61743098
/
sagecal
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

Merge pull request #41 from nlesc-dirac/sprint 

Sprint changes
This commit is contained in:
Faruk D 2018-02-28 15:49:20 +01:00 committed by GitHub
parent 22f89ae6d5 1088583a21
commit b4f101e0d6
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
23 changed files with 238 additions and 189 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

81
README.md
View File

@ -1,5 +1,17 @@
# SAGECAL
## Features
- Levenberg-Marquardt, LBFGS, Riemannian Trust Region, Nesterov's accelerated gradient descent algorithms
- GPU acceleration using CUDA
- Fast and accurate interferometric calibration
- Gaussian and Student's t noise models
- Shapelet source models
- CASA MS data format supported
- Distributed calibration using MPI - consensus optimization with data multiplexing
- Tools to build sky models and restore sky models to images
- Adaptive update of ADMM penalty (Barzilai-Borwein a.k.a. Spectral method)
Read INSTALL for installation. This file gives a brief guide to use SAGECal.
Warning: this file may be obsolete. use sagecal -h to see up-to-date options.
@ -7,76 +19,84 @@ Warning: this file may be obsolete. use sagecal -h to see up-to-date options.
## Step by Step Introduction:
###1a)Calibrate data in the standard way using BBS/CASA or anything else.
Use NDPP to average the data in your MS to a few channels (also average in time to about 10sec). Also flag any spikes in the data.
###1b)For subtraction of the ATeam from raw data (CasA,CygA,...), no initial calibration is necessary. Just run sagecal on raw data, but it is better to scale the sky model to match the apparent flux of the sources that are being subtracted.
### 1) Input Data
Input to sagecal must be in CASA MS format, make sure to create a column in the MS to write output data as well. The data can be in raw or averaged form, also initial calibration using other software can be also applied.
##2) Sky Model:
###2a)Make an image of your MS (using ExCon/casapy).
### 2) Sky Model:
#### 2a) Make an image of your MS (using ExCon/casapy).
Use Duchamp to create a mask for the image. Use buildsky to create a sky model. (see the README file on top level directory). Also create a proper cluster file.
Special options to buildsky: -o 1 (NOTE: not -o 2)
Alternatively, create these files by hand according to the following formats.
###2b)Cluster file format:
#### 2b) Cluster file format:
cluster_id chunk_size source1 source2 ...
e.g.
'''
```
0 1 P0C1 P0C2
2 3 P11C2 P11C1 P13C1
'''
```
Note: putting -ve values for cluster_id will not subtract them from data.
chunk_size: find hybrid solutions during one solve run. Eg. if -t 120 is used
to select 120 timeslots, cluster 0 will find a solution using the full 120 timeslots while cluster 2 will solve for every 120/3=40 timeslots.
###2c)Sky model format:
\#name h m s d m s I Q U V spectral_index RM extent_X(rad) extent_Y(rad) pos_angle(rad) freq0
#### 2c) Sky model format:
```
#name h m s d m s I Q U V spectral_index RM extent_X(rad) extent_Y(rad) pos_angle(rad) freq0
```
or
\#name h m s d m s I Q U V spectral_index1 spectral_index2 spectral_index3 RM extent_X(rad) extent_Y(rad) pos_angle(rad) freq0
```
#name h m s d m s I Q U V spectral_index1 spectral_index2 spectral_index3 RM extent_X(rad) extent_Y(rad) pos_angle(rad) freq0
```
e.g.
e.g.:
```
P1C1 0 12 42.996 85 43 21.514 0.030498 0 0 0 -5.713060 0 0 0 0 115039062.0
P5C1 1 18 5.864 85 58 39.755 0.041839 0 0 0 -6.672879 0 0 0 0 115039062.0
\#A Gaussian mjor,minor 0.1375,0.0917 deg diameter -> radius(rad), PA 43.4772 deg (-> rad)
\#Position Angle: "West from North (counter-clockwise)" (0 deg = North, 90 deg = West).
\#Note: PyBDSM and BBS use "North from East (counter-clockwise)" (0 deg = East, 90 deg = North).
#A Gaussian mjor,minor 0.1375,0.0917 deg diameter -> radius(rad), PA 43.4772 deg (-> rad)
#Position Angle: "West from North (counter-clockwise)" (0 deg = North, 90 deg = West).
#Note: PyBDSM and BBS use "North from East (counter-clockwise)" (0 deg = East, 90 deg = North).
G0 5 34 31.75 22 00 52.86 100 0 0 0 0.00 0 0.0012 0.0008 -2.329615801 130.0e6
\#A Disk radius=0.041 deg
#A Disk radius=0.041 deg
D01 23 23 25.67 58 48 58 80 0 0 0 0 0 0.000715 0.000715 0 130e6
\#A Ring radius=0.031 deg
#A Ring radius=0.031 deg
R01 23 23 25.416 58 48 57 70 0 0 0 0 0 0.00052 0.00052 0 130e6
\#A shapelet ('S3C61MD.fits.modes' file must be in the current directory)
#A shapelet ('S3C61MD.fits.modes' file must be in the current directory)
S3C61MD 2 22 49.796414 86 18 55.913266 0.135 0 0 0 -6.6 0 1 1 0.0 115000000.0
```
Note: Comments starting with a '#' are allowed for both sky model and cluster files.
Note: 3rd order spectral indices are also supported, use -F 1 option in sagecal.
Note: Spectral indices use natural logarithm, exp(ln(I0) + p1*ln(f/f0) + p2*ln(f/f0)^2 + ..) so if you have a model with common logarithms like 10^(log(J0) + q1*log(f/f0) + q2*log(f/f0)^2 + ..) then, conversion is
Note: Spectral indices use natural logarithm, ```exp(ln(I0) + p1 * ln(f/f0) + p2 * ln(f/f0)^2 + ..)``` so if you have a model with common logarithms like ```10^(log(J0) + q1*log(f/f0) + q2*log(f/f0)^2 + ..)``` then, conversion is
```
ln(I0)+p1*ln(f/f0)+p2*ln(f/f0)^2+... = ln(10)*(log(J0)+q1*log(f/f0)+q2*log(f/f0))^2)+...)
=ln(10)*(ln(J0)/ln(10)+q1*ln(f/f0)/ln(10)+q2*ln(f/f0)^2/ln(10)^2+...)
```
so
```
I0=J0
p1=q1
p2=q2/ln(10)
p3=q3/ln(10)^2
...
```
##3)Run sagecal
### 3) Run sagecal
Optionally: Make sure your machine has (1/2 working NVIDIA GPU cards or Intel Xeon Phi MICs) to use sagecal.
Recommended usage: (with GPUs)
```
sagecal -d my_data.MS -s my_skymodel -c my_clustering -n no.of.threads -t 60 -p my_solutions -e 3 -g 2 -l 10 -m 7 -w 1 -b 1
```
Use your solution interval (-t 60) so that its big enough to get a decent solution and not too big to make the parameters vary too much. (about 20 minutes per solution is reasonable).
@ -89,20 +109,24 @@ With -a 1 and -p 'solutions_file', simulation is done with the sky model corrupt
With -a 1 and -p 'solutions_file' and -z 'ignore_file', simulation is done with the solutions in the 'solutions_file', but ignoring the cluster ids in the 'ignore_file'.
Eg. If you need to ignore cluster ids '-1', '10', '999', create a text file :
```
-1
10
999
```
and use it as the 'ignore_file'.
##4)Distributed calibration
### 4) Distributed calibration
Use mpirun to run sagecal-mpi, example:
```
mpirun -np 11 -hostfile ./machines --map-by node --cpus-per-proc 8
--mca yield_when_idle 1 -mca orte_tmpdir_base /scratch/users/sarod
/full/path/to/sagecal-mpi -f 'MS*pattern' -A 30 -P 2 -r 5
-s sky.txt -c cluster.txt -n 16 -t 1 -e 3 -g 2 -l 10 -m 7 -x 10 -F 1 -j 5
```
Specific options :
-np 11 : 11 processes : starts 10 slaves + 1 master
@ -119,7 +143,7 @@ Note: the number of slaves (-np option) can be lower than the number of MS calib
The rest of the options are similar to sagecal.
##5)Solution format
### 5) Solution format
All SAGECal solutions are stored as text files. Lines starting with '#' are comments.
The first non-comment line includes some general information, i.e.
freq(MHz) bandwidth(MHz) time_interval(min) stations clusters effective_clusters
@ -128,7 +152,12 @@ The remaining lines contain solutions for each cluster as a single column, the f
Let's say there are K effective clusters and N directions. Then there will be K+1 columns, the first column will start from 0 and increase to 8N-1,
which can be used to count the row number. It will keep repeating this, for each time interval.
The rows 0 to 7 belong to the solutions for the 1st station. The rows 8 to 15 for the 2nd station and so on.
Each 8 rows of any given column represent the 8 values of a 2x2 Jones matrix. Lets say these are S0,S1,S2,S3,S4,S5,S6 and S7. Then the Jones matrix is [S0+j*S1, S4+j*S5; S2+j*S3, S6+j*S7] (the ';' denotes the 1st row of the 2x2 matrix).
Each 8 rows of any given column represent the 8 values of a 2x2 Jones matrix. Lets say these are ```S0,S1,S2,S3,S4,S5,S6``` and ```S7```. Then the Jones matrix is ```[S0+j*S1, S4+j*S5; S2+j*S3, S6+j*S7]``` (the ';' denotes the 1st row of the 2x2 matrix).
When a luster has a chunk size > 1, there will be more than 1 solution per given time interval.
So for this cluster, there will be more than 1 column in the solution file, the exact number of columns being equal to the chunk size.
### Additional Info
See [LOFAR Cookbook Chapter](https://support.astron.nl/LOFARImagingCookbook/sagecal.html).

2
src/MPI/sagecal_master.cpp
View File

@ -403,7 +403,7 @@ cout<<"Reference frequency (MHz)="<<iodata.freq0*1.0e-6<<endl;
FILE *dfp;
if ((dfp=fopen("debug.m","w+"))==0) {
fprintf(stderr,"%s: %d: no file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
fprintf(dfp,"%% timeslot admmIter dir1 dir2 (reversed) ...; (each row=freq)\n");
fprintf(dfp,"rhoK=[\n");

2
src/MPI/sagecal_slave.cpp
View File

@ -262,7 +262,7 @@ cout<<"Slave "<<myrank<<" has nothing to do"<<endl;
FILE *sfq;
if ((sfq=fopen(initsolfile,"r"))==0) {
fprintf(stderr,"%s: %d: no solution file present\n",__FILE__,__LINE__);
return 1;
exit(1);
}
/* remember to skip first 3 lines from solution file */
char chr;

94
src/buildsky/buildmultisky.c
View File

@ -350,27 +350,27 @@ read_fits_file_f(const char *fitsdir, const char *maskfile, GHashTable **pixtabl
ncoord=g_list_length(val->pix);
if ((pixelc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((worldc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((phic=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((thetac=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((statc=(int*)calloc((size_t)ncoord,sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
kk=0;
@ -463,27 +463,27 @@ read_fits_file_f(const char *fitsdir, const char *maskfile, GHashTable **pixtabl
ncoord=1; /* memory for WCS, 1 pixel */
if ((pixelc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((worldc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((phic=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((thetac=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((statc=(int*)calloc((size_t)ncoord,sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
kk=0;
pixelc[kk+0]=(double)1.0;
@ -829,15 +829,15 @@ add_guard_pixels_f(GList *pixlist, int Nf, double threshold, hpixelf **parr, int
/* insert this pixel to hash table */
if ((xykey = (xyhash*)malloc(sizeof(xyhash)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((ykey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
@ -922,23 +922,23 @@ add_guard_pixels_f(GList *pixlist, int Nf, double threshold, hpixelf **parr, int
deltav=((*mp-*lp)/(double)origx);
if ((lv= (double*)malloc(sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*lv=*lp-deltav;
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*xkey=xarr[0];
g_hash_table_insert(xtol,(gpointer)xkey,(gpointer)lv);
if ((lv= (double*)malloc(sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*lv=*mp+deltav;
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*xkey=xarr[origx+1];
g_hash_table_insert(xtol,(gpointer)xkey,(gpointer)lv);
@ -956,7 +956,7 @@ add_guard_pixels_f(GList *pixlist, int Nf, double threshold, hpixelf **parr, int
/* SORT */
if ((yarr= (uint32_t*)malloc(sizeof(uint32_t)*(origy+2)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
pixval0=pixval;
ci=1;
@ -983,23 +983,23 @@ add_guard_pixels_f(GList *pixlist, int Nf, double threshold, hpixelf **parr, int
deltav=((*mp-*lp)/(double)origy);
if ((lv= (double*)malloc(sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*lv=*lp-deltav;
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*xkey=yarr[0];
g_hash_table_insert(ytom,(gpointer)xkey,(gpointer)lv);
if ((lv= (double*)malloc(sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*lv=*mp+deltav;
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*xkey=yarr[origy+1];
g_hash_table_insert(ytom,(gpointer)xkey,(gpointer)lv);
@ -1016,7 +1016,7 @@ add_guard_pixels_f(GList *pixlist, int Nf, double threshold, hpixelf **parr, int
*n=g_hash_table_size(xtol)*g_hash_table_size(ytom);
if ((*parr=(hpixelf*)malloc(sizeof(hpixelf)*(*n)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
minflux=1e6;
ck=0;
@ -1039,7 +1039,7 @@ add_guard_pixels_f(GList *pixlist, int Nf, double threshold, hpixelf **parr, int
(*parr)[ck].Nf=Nf;
if (((*parr)[ck].sI = (double*)malloc(sizeof(double)*Nf))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
memcpy((*parr)[ck].sI,ppix->sI,sizeof(double)*Nf);
for (ci=0; ci<Nf; ci++) {
@ -1056,15 +1056,15 @@ add_guard_pixels_f(GList *pixlist, int Nf, double threshold, hpixelf **parr, int
/* now the guard pixels */
if ((xykey = (xyhash*)malloc(sizeof(xyhash)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((ykey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
minflux*=threshold;
@ -1088,7 +1088,7 @@ add_guard_pixels_f(GList *pixlist, int Nf, double threshold, hpixelf **parr, int
(*parr)[ck].Nf=Nf;
if (((*parr)[ck].sI= (double*)malloc(sizeof(double)*(Nf)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
for (cl=0; cl<Nf; ++cl) {
(*parr)[ck].sI[cl]=minflux;
@ -1363,7 +1363,7 @@ write_world_coords_f(const char *imgfile, GHashTable *pixtable, double minpix, i
/* read FITS header */
if ((status = fits_hdr2str(fbuff.fptr, 1, NULL, 0, &header, &ncard, &status))) {
fits_report_error(stderr, status);
return 1;
exit(1);
}
/* try to Parse the primary header of the FITS file. */
@ -1382,54 +1382,54 @@ write_world_coords_f(const char *imgfile, GHashTable *pixtable, double minpix, i
if ((status = wcsset(fbuff.wcs))) {
fprintf(stderr, "wcsset ERROR %d:\n", status);
return 1;
//return 1; not something serious, just a warning is fine
}
ncoord=1;
if ((pixell=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((pixelm=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgl=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgm=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((ra_c=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((dec_c=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((statc=(int*)calloc((size_t)ncoord,sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
outf=fopen(textfile,"w+");
if(!outf) {
fprintf(stderr,"%s: %d: unable to open file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
regf=fopen(regionfile,"w+");
if(!regf) {
fprintf(stderr,"%s: %d: unable to open file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
flxf=fopen(fluxfile,"w+");
if(!flxf) {
fprintf(stderr,"%s: %d: unable to open file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if (!outformat) {
fprintf(outf,"# (Name, Type, Ra, Dec, I, Q, U, V, ReferenceFrequency='%12.7lg', SpectralIndex='[0.0,0.0,0.0]', MajorAxis, MinorAxis, Orientation) = format\n",ref_freq);
@ -1802,12 +1802,12 @@ filter_pixels_f(GHashTable *pixtable, double wcutoff) {
/* A: size 2x2 */
if ((A=(double*)malloc(sizeof(double)*(size_t)N*N))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
/* W: size 2x1 */
if ((W=(double*)malloc(sizeof(double)*(size_t)N))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
printf("################## probable ignore list (%lf) ##############\n",wcutoff);
@ -1825,15 +1825,15 @@ filter_pixels_f(GHashTable *pixtable, double wcutoff) {
Npix=g_list_length(val->pix);
if ((xpix=(double*)malloc(sizeof(double)*(size_t)Npix))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((ypix=(double*)malloc(sizeof(double)*(size_t)Npix))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((sI=(double*)malloc(sizeof(double)*(size_t)Npix))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
ci=0;
xmean=ymean=0.0;

76
src/buildsky/buildsky.c
View File

@ -432,27 +432,27 @@ read_fits_file(const char *imgfile, const char *maskfile, GHashTable **pixtable,
ncoord=g_list_length(val->pix);
if ((pixelc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((worldc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((phic=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((thetac=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((statc=(int*)calloc((size_t)ncoord,sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
kk=0;
@ -629,7 +629,7 @@ write_world_coords(const char *imgfile, GHashTable *pixtable, double minpix, dou
/* read FITS header */
if ((status = fits_hdr2str(fbuff.fptr, 1, NULL, 0, &header, &ncard, &status))) {
fits_report_error(stderr, status);
return 1;
exit(1);
}
/* try to Parse the primary header of the FITS file. */
@ -648,49 +648,49 @@ write_world_coords(const char *imgfile, GHashTable *pixtable, double minpix, dou
if ((status = wcsset(fbuff.wcs))) {
fprintf(stderr, "wcsset ERROR %d:\n", status);
return 1;
//return 1; just a warning is enough here
}
ncoord=1;
if ((pixell=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((pixelm=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgl=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgm=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((ra_c=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((dec_c=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((statc=(int*)calloc((size_t)ncoord,sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
outf=fopen(textfile,"w+");
if(!outf) {
fprintf(stderr,"%s: %d: unable to open file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
regf=fopen(regionfile,"w+");
if(!regf) {
fprintf(stderr,"%s: %d: unable to open file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if (!outformat) {
fprintf(outf,"# (Name, Type, Ra, Dec, I, Q, U, V, ReferenceFrequency='60e6', SpectralIndex='[0.0,0.0,0.0]', MajorAxis, MinorAxis, Orientation) = format\n");
@ -1011,15 +1011,15 @@ add_guard_pixels(GList *pixlist, double threshold, hpixel **parr, int *n) {
/* insert this pixel to hash table */
if ((xykey = (xyhash*)malloc(sizeof(xyhash)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((ykey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
@ -1104,23 +1104,23 @@ add_guard_pixels(GList *pixlist, double threshold, hpixel **parr, int *n) {
deltav=((*mp-*lp)/(double)origx);
if ((lv= (double*)malloc(sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*lv=*lp-deltav;
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*xkey=xarr[0];
g_hash_table_insert(xtol,(gpointer)xkey,(gpointer)lv);
if ((lv= (double*)malloc(sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*lv=*mp+deltav;
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*xkey=xarr[origx+1];
g_hash_table_insert(xtol,(gpointer)xkey,(gpointer)lv);
@ -1138,7 +1138,7 @@ add_guard_pixels(GList *pixlist, double threshold, hpixel **parr, int *n) {
/* SORT */
if ((yarr= (uint32_t*)malloc(sizeof(uint32_t)*(origy+2)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
pixval0=pixval;
ci=1;
@ -1165,23 +1165,23 @@ add_guard_pixels(GList *pixlist, double threshold, hpixel **parr, int *n) {
deltav=((*mp-*lp)/(double)origy);
if ((lv= (double*)malloc(sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*lv=*lp-deltav;
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*xkey=yarr[0];
g_hash_table_insert(ytom,(gpointer)xkey,(gpointer)lv);
if ((lv= (double*)malloc(sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*lv=*mp+deltav;
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*xkey=yarr[origy+1];
g_hash_table_insert(ytom,(gpointer)xkey,(gpointer)lv);
@ -1198,7 +1198,7 @@ add_guard_pixels(GList *pixlist, double threshold, hpixel **parr, int *n) {
*n=g_hash_table_size(xtol)*g_hash_table_size(ytom);
if ((*parr= (hpixel *)malloc(sizeof(hpixel)*(*n)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
minflux=1e6;
ck=0;
@ -1227,15 +1227,15 @@ add_guard_pixels(GList *pixlist, double threshold, hpixel **parr, int *n) {
/* now the guard pixels */
if ((xykey = (xyhash*)malloc(sizeof(xyhash)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((xkey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((ykey = (uint32_t*)malloc(sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
minflux*=threshold;
@ -1452,12 +1452,12 @@ filter_pixels(GHashTable *pixtable, double wcutoff) {
/* A: size 2x2 */
if ((A=(double*)malloc(sizeof(double)*(size_t)N*N))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
/* W: size 2x1 */
if ((W=(double*)malloc(sizeof(double)*(size_t)N))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
printf("################## probable ignore list (%lf) ##############\n",wcutoff);
@ -1475,15 +1475,15 @@ filter_pixels(GHashTable *pixtable, double wcutoff) {
Npix=g_list_length(val->pix);
if ((xpix=(double*)malloc(sizeof(double)*(size_t)Npix))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((ypix=(double*)malloc(sizeof(double)*(size_t)Npix))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((sI=(double*)malloc(sizeof(double)*(size_t)Npix))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
ci=0;
xmean=ymean=0.0;

15
src/buildsky/cluster.c
View File

@ -2654,7 +2654,7 @@ number of clusters is larger than the number of elements being clustered,
int* mapping = NULL;
double** cdata;
int** cmask;
int* counts;
int* counts=0;
if (nelements < nclusters)
{ *ifound = 0;
@ -2691,6 +2691,7 @@ number of clusters is larger than the number of elements being clustered,
else ok = makedatamask(ndata, nclusters, &cdata, &cmask);
if(!ok)
{ free(counts);
counts=0;
if(npass>1)
{ free(tclusterid);
free(mapping);
@ -2720,8 +2721,7 @@ number of clusters is larger than the number of elements being clustered,
if (transpose==0) freedatamask(nclusters, cdata, cmask);
else freedatamask(ndata, cdata, cmask);
free(counts);
if (counts) free(counts);
}
/* *********************************************************************** */
@ -2990,6 +2990,7 @@ when microarrays are being clustered.
if (i < n) /* break condition encountered */
{ j = i;
for (i = 1; i < j; i++) free(matrix[i]);
free(matrix);
return NULL;
}
@ -3526,9 +3527,13 @@ If a memory error occurs, pslcluster returns NULL.
free(vector);
free(index);
result = realloc(result, nnodes*sizeof(Node));
return result;
//result = realloc(result, nnodes*sizeof(Node));
Node* result1 = malloc(nnodes*sizeof(Node)); /* nnodes = nlements -1 */
memcpy(result1,result,nnodes*sizeof(Node));
free(result);
return result1;
}
/* ******************************************************************** */

4
src/buildsky/hull.c
View File

@ -204,7 +204,7 @@ construct_boundary(pixellist *pixset)
pixset->Nh=count(&sl)+count(&su)-1;
if ((pixset->hull= (hpoint*)malloc((size_t)pixset->Nh*sizeof(hpoint)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
i=0;
@ -343,7 +343,7 @@ construct_boundary_f(pixellistf *pixset)
pixset->Nh=count(&sl)+count(&su)-1;
if ((pixset->hull= (hpoint*)malloc((size_t)pixset->Nh*sizeof(hpoint)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
i=0;

1
src/buildsky/main.c
View File

@ -310,6 +310,7 @@ int main(int argc, char **argv) {
}
/* if filter is on, quit now */
printf("quitting. re-run without filter\n");
free(unistr);
return 0;
}
if (!multifits) {

32
src/buildsky/scluster.c
View File

@ -114,24 +114,24 @@ fit_two_to_one(GList *pixlist, double bmaj, double bmin, double bpa, double *lva
if ((pixell=(double*)calloc((size_t)npix,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((pixelm=(double*)calloc((size_t)npix,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((p=(double*)calloc((size_t)m,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((x=(double*)calloc((size_t)n,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((x1=(double*)calloc((size_t)n,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
ci=0;
for(pli=pixlist; pli!=NULL; pli=g_list_next(pli)) {
@ -267,24 +267,24 @@ fit_two_to_one_f(GList *pixlist, int Nf, double *freqs, double *bmaj, double *bm
if ((pixell=(double*)calloc((size_t)npix,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((pixelm=(double*)calloc((size_t)npix,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((p=(double*)calloc((size_t)m,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((x=(double*)calloc((size_t)n,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((x1=(double*)calloc((size_t)n,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
ci=0;
for(pli=pixlist; pli!=NULL; pli=g_list_next(pli)) {
@ -389,12 +389,12 @@ cluster_sources(double r, GList *inlist, GList *pixlist, double bmaj, double bmi
for (li=inlist; li!=NULL; li=g_list_next(li)) {
if((qnode= (pqsrc*)malloc(sizeof(pqsrc)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
qnode->rd=0.0;
if((qnode->src= (extsrc*)malloc(sizeof(extsrc)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
memcpy((void*)qnode->src,(void*)li->data,sizeof(extsrc)); /* copy source data */
@ -541,12 +541,12 @@ cluster_sources_f(double r, GList *inlist, GList *pixlist, int Nf, double *freqs
for (li=inlist; li!=NULL; li=g_list_next(li)) {
if((qnode= (pqsrcf*)malloc(sizeof(pqsrcf)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
qnode->rd=0.0;
if((qnode->src=(extsrcf*)malloc(sizeof(extsrcf)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
memcpy((void*)qnode->src,(void*)li->data,sizeof(extsrcf)); /* copy source data */
@ -887,12 +887,12 @@ cluster_sky(const char *imgfile, GList *skylist, int ncluster) {
cfilep=fopen(clusterfile,"w+");
if(!cfilep) {
fprintf(stderr,"%s: %d: unable to open file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
annfp=fopen(clusterann,"w+");
if(!annfp) {
fprintf(stderr,"%s: %d: unable to open file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
/* print header info in the files */

4
src/lib/Dirac/clmfit_nocuda.c
View File

@ -878,7 +878,7 @@ mlm_der_single(
printf("Singular matrix info=%d\n",status);
#endif
}
} else {
} else if (solve_axb==2) {
/* SVD solver *********************************/
/* U S VT = A */
status=my_dgesvd('A','A',M,M,JkTJk,M,Sd,Ud,M,VTd,M,WORK,lwork);
@ -940,7 +940,7 @@ mlm_der_single(
printf("Singular matrix info=%d\n",status);
#endif
}
} else {
} else if (solve_axb==2) {
/* SVD solver *********************************/
/* dhatk <= U^T jacTed */
my_dgemv('T',M,M,1.0,Ud,M,JkTe,1,0.0,dhatk,1);

4
src/lib/Dirac/consensus_poly.c
View File

@ -190,7 +190,7 @@ setup_polynomials(double *B, int Npoly, int Nf, double *freqs, double freq0, int
int
find_prod_inverse(double *B, double *Bi, int Npoly, int Nf, double *fratio) {
int ci,status,lwork=0;
int ci,status,lwork=1;
double w[1],*WORK,*U,*S,*VT;
/* set Bi to zero */
memset(Bi,0,sizeof(double)*Npoly*Npoly);
@ -301,7 +301,7 @@ sum_inv_threadfn(void *data) {
thread_data_prod_inv_t *t=(thread_data_prod_inv_t*)data;
double w[1],*WORK,*U,*S,*VT;
int k,ci,status,lwork=0;
int k,ci,status,lwork=1;
int Np2=t->Npoly*t->Npoly;
/* allocate memory for the SVD here */
if ((U=(double*)calloc((size_t)Np2,sizeof(double)))==0) {

2
src/lib/Dirac/lbfgs.c
View File

@ -225,7 +225,7 @@ mult_hessian(int m, double *pk, double *gk, double *s, double *y, double *rho, i
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
exit(1);
}
if ((idx=(int*)calloc((size_t)M,sizeof(double)))==0) {
if ((idx=(int*)calloc((size_t)M,sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
exit(1);
}

2
src/lib/Dirac/lbfgs_nocuda.c
View File

@ -318,7 +318,7 @@ mult_hessian(int m, double *pk, double *gk, double *s, double *y, double *rho, i
#endif
exit(1);
}
if ((idx=(int*)calloc((size_t)M,sizeof(double)))==0) {
if ((idx=(int*)calloc((size_t)M,sizeof(int)))==0) {
#ifndef USE_MIC
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
#endif

6
src/lib/Dirac/mdl.c
View File

@ -63,7 +63,7 @@ minimum_description_length(int N, int M, int F, double *J, double *rho, double *
int m,q;
if ((dfp=fopen("debug.m","w+"))==0) {
fprintf(stderr,"%s: %d: no file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
for(cm=0; cm<F; cm++) {
@ -234,7 +234,9 @@ minimum_description_length(int N, int M, int F, double *J, double *rho, double *
/* MDL = F/2 log(RSS/F) + Npoly/2 log(F) */
mdl[idx++]=(0.5*(double)F)*log(RSS/(double)(F))+0.5*(double)(Npoly)*log((double)F);
#ifdef DEBUG
fclose(dfp);
#endif
free(Z);
free(z);

2
src/lib/Dirac/robust_lbfgs_nocuda.c
View File

@ -165,7 +165,7 @@ mult_hessian(int m, double *pk, double *gk, double *s, double *y, double *rho, i
#endif
exit(1);
}
if ((idx=(int*)calloc((size_t)M,sizeof(double)))==0) {
if ((idx=(int*)calloc((size_t)M,sizeof(int)))==0) {
#ifndef USE_MIC
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
#endif

2
src/lib/Dirac/rtr_solve_robust.c
View File

@ -1472,7 +1472,7 @@ itrr(int N,complex double *x,complex double *eta, complex double *Heta, double *
sigma=delta_0;
}
/* radius update */
double beta_1,beta_2,beta_i;
double beta_1,beta_2,beta_i=0.0;
beta_1=0.0;
beta_2=0.0;
if (m<MK) {

2
src/lib/Dirac/rtr_solve_robust_admm.c
View File

@ -1456,7 +1456,7 @@ itrr(int N,complex double *x,complex double *eta, complex double *Heta, double *
sigma=delta_0;
}
/* radius update */
double beta_1,beta_2,beta_i;
double beta_1,beta_2,beta_i=0.0;
beta_1=0.0;
beta_2=0.0;
if (m<MK) {

2
src/lib/Dirac/rtr_solve_robust_cuda.c
View File

@ -361,7 +361,7 @@ itrr(int ThreadsPerBlock, int BlocksPerGrid, int N, int M, cuFloatComplex *x, cu
sigma=delta_0;
}
/* radius update */
float beta_1,beta_2,beta_i;
float beta_1,beta_2,beta_i=0.0f;
beta_1=0.0f;
beta_2=0.0f;

2
src/lib/Dirac/rtr_solve_robust_cuda_admm.c
View File

@ -274,7 +274,7 @@ itrr(int ThreadsPerBlock, int BlocksPerGrid, int N, int M, cuFloatComplex *x, cu
sigma=delta_0;
}
/* radius update */
float beta_1,beta_2,beta_i;
float beta_1,beta_2,beta_i=0.0f;
beta_1=0.0f;
beta_2=0.0f;

5
src/lib/Radio/readsky.c
View File

@ -177,6 +177,7 @@ read_shapelet_modes(char *buff,int *n0,double *beta,double **modes) {
}
}
fclose(cfp);
free(input_modes);
return 0;
}
@ -685,8 +686,10 @@ printf("Nc=%d\n",Nc);
memset(buf,0,buff_len);
c=read_next_string(&buf,&buff_len,cfp);
if (c!=1) {
/* first column is solution number (int) 1..8N */
/* first column is solution number (int) 0..8N-1 */
sscanf(buf,"%d",&cn);
/* also do a sanity check */
if ( cn<0 || cn>Nc ) cn=0;
}
#ifdef DEBUG
printf("%d ",cn);

1
src/lib/Radio/residual.c
View File

@ -1714,6 +1714,7 @@ predict_visibilities_multifreq_withsol(double *u,double *v,double *w,double *p,d
pthread_attr_destroy(&attr);
if (pinv) free(pinv);
free(th_array);
free(threaddata);

31
src/restore/readsky.c
View File

@ -173,6 +173,7 @@ read_shapelet_modes_img(char *buff,int *n0,double *beta,double **modes, double *
}
*eX=l_a; *eY=l_b; *eP=l_theta; *ltf=linear_tf;
fclose(cfp);
free(input_modes);
return 0;
}
@ -520,7 +521,7 @@ read_sky_model_withgain(const char *slistname, glist *slist, int format,const ch
/* new cluster found */
if ((clus=(clust_t*)malloc(sizeof(clust_t)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
sscanf(buf,"%d",&clus->id);
clus->slist=NULL;
@ -549,11 +550,11 @@ read_sky_model_withgain(const char *slistname, glist *slist, int format,const ch
/* source found for this cluster */
if ((sclus= (clust_n*)malloc(sizeof(clust_n)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((sclus->name=(char*)malloc((size_t)(strlen(buf)+1)*sizeof(char)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
strcpy(sclus->name,buf);
#ifdef DEBUG
@ -624,12 +625,12 @@ read_sky_model_withgain(const char *slistname, glist *slist, int format,const ch
sclus=ln->data;
if ((hkey=(char*)malloc((size_t)(strlen(sclus->name)+1)*sizeof(char)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
strcpy(hkey,sclus->name);
if ((cptr=(clust_t**)malloc((size_t)(1)*sizeof(clust_t*)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
*cptr=clus;
g_hash_table_insert(stable,(gpointer)hkey,(gpointer)cptr);
@ -679,38 +680,38 @@ read_sky_model_withgain(const char *slistname, glist *slist, int format,const ch
/* memory to store solutions */
if ((JHJ=(complex double**)calloc((size_t)(4),sizeof(complex double*)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((JsHJ=(complex double**)calloc((size_t)(4),sizeof(complex double*)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
for (ci=0; ci<4; ci++) {
if ((JHJ[ci]=(complex double*)calloc((size_t)(Nc),sizeof(complex double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((JsHJ[ci]=(complex double*)calloc((size_t)(Nc),sizeof(complex double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
}
if ((Jf1=(double**)calloc((size_t)(8),sizeof(double*)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((Jf2=(double**)calloc((size_t)(8),sizeof(double*)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
for (ci=0; ci<8; ci++) {
if ((Jf1[ci]=(double*)calloc((size_t)(Nc),sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((Jf2[ci]=(double*)calloc((size_t)(Nc),sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
}
@ -720,7 +721,7 @@ read_sky_model_withgain(const char *slistname, glist *slist, int format,const ch
*/
if ((cfp=fopen(solfile,"r"))==0) {
fprintf(stderr,"%s: %d: no file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
/* for new versions of sagecal, need to skip first 3 lines */
for (ci=0; ci<3; ci++) {
@ -843,7 +844,7 @@ read_sky_model_withgain(const char *slistname, glist *slist, int format,const ch
/* read in the source list, parse each source info */
if ((cfp=fopen(slistname,"r"))==0) {
fprintf(stderr,"%s: %d: no file\n",__FILE__,__LINE__);
return 1;
exit(1);
}
c=skip_lines(cfp);

55
src/restore/restore.c
View File

@ -55,9 +55,10 @@ zero_image(long totalrows, long offset, long firstrow, long nrows,
charp= (float *) fits_iter_get_array(&cols[0]);
}
/* NOTE: 1st element of array is the null pixel value! */
/* Loop from 1 to nrows, not 0 to nrows - 1. */
for (ii = 1; ii <= nrows; ii++) {
charp[ii]=0.0;
charp[ii]=0.0;
}
return 0;
@ -88,32 +89,32 @@ calculate_contribution1(struct wcsprm *wcs, double ll, double mm, sinfo *ss, dou
ncoord=1;
if ((pixelra=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((pixeldec=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgphi=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgtheta=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgl=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgm=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((statc=(int*)calloc((size_t)ncoord,sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
pixelra[0]=(double)(ss->ra)*180.0/M_PI;
pixeldec[0]=(double)(ss->dec)*180.0/M_PI;
@ -433,27 +434,27 @@ read_fits_file_restore(const char *filename, glist *slist, double bmaj,double bm
ncoord=1;
if ((pixelc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((worldc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((phic=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((thetac=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((statc=(int*)calloc((size_t)ncoord,sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
/* if format==FORMAT_LM, convert l,m of source to ra,dec */
if (format==FORMAT_LM) {
@ -621,34 +622,34 @@ read_fits_file_restore(const char *filename, glist *slist, double bmaj,double bm
if ((myarr=(double*)calloc((size_t)totalpix,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
/* allocate memory for pixel/world coordinate arrays */
ncoord=new_naxis[0]*new_naxis[1]*1*1; /* consider only one plane fron freq, and stokes axes because RA,Dec will not change */
if ((pixelc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((imgc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((worldc=(double*)calloc((size_t)ncoord*4,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((phic=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((thetac=(double*)calloc((size_t)ncoord,sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((statc=(int*)calloc((size_t)ncoord,sizeof(int)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
/* fill up the pixel coordinate array */
@ -763,11 +764,11 @@ read_fits_file_restore(const char *filename, glist *slist, double bmaj,double bm
M=exs->n0*exs->n0-1; /* always less than */
if ((x=(double*)calloc((size_t)new_naxis[0],sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
if ((y=(double*)calloc((size_t)new_naxis[1],sizeof(double)))==0) {
fprintf(stderr,"%s: %d: no free memory\n",__FILE__,__LINE__);
return 1;
exit(1);
}
Nx=new_naxis[0];
Ny=new_naxis[1];
@ -1038,6 +1039,12 @@ int main(int argc, char **argv) {
free(clusterfile);
glist_delete(&slist);
free(ignfile);
} else {
if (ffile) free(ffile);
if (slistname) free(slistname);
if (solfile) free(solfile);
if (clusterfile) free(clusterfile);
if (ignfile) free(ignfile);
}
return 0;
}