|
| 1 | +/* Angle gather at fixed lateral position in 3-D */ |
| 2 | +/* |
| 3 | + Copyright (C) 2012 University of Texas at Austin |
| 4 | +
|
| 5 | + This program is free software; you can redistribute it and/or modify |
| 6 | + it under the terms of the GNU General Public License as published by |
| 7 | + the Free Software Foundation; either version 2 of the License, or |
| 8 | + (at your option) any later version. |
| 9 | +
|
| 10 | + This program is distributed in the hope that it will be useful, |
| 11 | + but WITHOUT ANY WARRANTY; without even the implied warranty of |
| 12 | + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
| 13 | + GNU General Public License for more details. |
| 14 | +
|
| 15 | + You should have received a copy of the GNU General Public License |
| 16 | + along with this program; if not, write to the Free Software |
| 17 | + Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA |
| 18 | +*/ |
| 19 | + |
| 20 | +#include <rsf.h> |
| 21 | + |
| 22 | +#ifndef _cram_gather3_h |
| 23 | + |
| 24 | +#include "cram_point3.h" |
| 25 | +/*^*/ |
| 26 | + |
| 27 | +typedef struct CRAMGather3 *sf_cram_gather3; |
| 28 | +/* abstract data type */ |
| 29 | +/*^*/ |
| 30 | + |
| 31 | +#endif |
| 32 | + |
| 33 | +struct CRAMGather3 { |
| 34 | + int nb, na, nz, noa, nda, oda; |
| 35 | + float b0, a0, db, da; |
| 36 | + float *stack; |
| 37 | + float *energy; |
| 38 | + float *illum; |
| 39 | + bool outaz, inorm; |
| 40 | +}; |
| 41 | +/* concrete data type */ |
| 42 | + |
| 43 | +sf_cram_gather3 sf_cram_gather3_init (int nb, int na, int nz, |
| 44 | + float b0, float db, float a0, float da, |
| 45 | + float oamax, float damax, bool outaz, bool inorm) |
| 46 | +/*< Initialize object >*/ |
| 47 | +{ |
| 48 | + sf_cram_gather3 cram_gather = (sf_cram_gather3)sf_alloc (1, sizeof (struct CRAMGather3)); |
| 49 | + |
| 50 | + cram_gather->nb = nb; /* Inclination */ |
| 51 | + cram_gather->b0 = b0; |
| 52 | + cram_gather->db = db; |
| 53 | + cram_gather->na = na; /* Azimuth */ |
| 54 | + cram_gather->a0 = a0; |
| 55 | + cram_gather->da = da; |
| 56 | + cram_gather->nz = nz; /* Number of depth samples */ |
| 57 | + |
| 58 | + cram_gather->outaz = outaz; /* Whether to output azimuth dimension */ |
| 59 | + cram_gather->inorm = inorm; /* Whether to normalize gathers */ |
| 60 | + |
| 61 | + /* Maximum number of scattering angles in the output */ |
| 62 | + cram_gather->noa = (int)(oamax/db + 0.5) + 1; |
| 63 | + if (cram_gather->noa > nb) |
| 64 | + cram_gather->noa = nb; |
| 65 | + /* Maximum number of dip angles in the output */ |
| 66 | + cram_gather->nda = (int)(2*damax/(0.5*db) + 0.5) + 1; |
| 67 | + if (cram_gather->nda > 4*nb) |
| 68 | + cram_gather->nda = 4*nb; |
| 69 | + /* Index of the first dip angle in the output */ |
| 70 | + cram_gather->oda = (4*nb - cram_gather->nda)/2; |
| 71 | + if (cram_gather->oda < 0) |
| 72 | + cram_gather->oda = 0; |
| 73 | + |
| 74 | + /* Buffers for collapsing azimuth dimension */ |
| 75 | + cram_gather->stack = sf_floatalloc (4*nb); |
| 76 | + cram_gather->energy = sf_floatalloc (4*nb); |
| 77 | + cram_gather->illum = sf_floatalloc (4*nb); |
| 78 | + |
| 79 | + return cram_gather; |
| 80 | +} |
| 81 | + |
| 82 | +void sf_cram_gather3_close (sf_cram_gather3 cram_gather) |
| 83 | +/*< Destroy object >*/ |
| 84 | +{ |
| 85 | + free (cram_gather->stack); |
| 86 | + free (cram_gather->energy); |
| 87 | + free (cram_gather->illum); |
| 88 | + free (cram_gather); |
| 89 | +} |
| 90 | + |
| 91 | +void sf_cram_gather3_image_output (sf_cram_gather3 cram_gather, sf_cram_point3 cram_point, |
| 92 | + sf_file imfile, sf_file ilfile) |
| 93 | +/*< Write one spatial location from the image to disk >*/ |
| 94 | +{ |
| 95 | + float smp, hits; |
| 96 | + |
| 97 | + smp = sf_cram_point3_get_image (cram_point, &hits); |
| 98 | + sf_floatwrite (&smp, 1, imfile); |
| 99 | + if (ilfile) |
| 100 | + sf_floatwrite (&hits, 1, ilfile); |
| 101 | +} |
| 102 | + |
| 103 | +void sf_cram_gather3_oangle_output (sf_cram_gather3 cram_gather, sf_cram_point3 cram_point, |
| 104 | + sf_file oafile, sf_file sqfile, sf_file ilfile) |
| 105 | +/*< Add one spatial point to the scattering angle gather and write it to disk >*/ |
| 106 | +{ |
| 107 | + int ib, ia; |
| 108 | + float **oimage, **osqimg, **ohits; |
| 109 | + |
| 110 | + oimage = sf_cram_point3_get_oimage (cram_point, &osqimg, &ohits); |
| 111 | + |
| 112 | + if (cram_gather->inorm) { |
| 113 | + for (ia = 0; ia < cram_gather->na/4; ia++) { |
| 114 | + for (ib = 0; ib < 2*cram_gather->nb; ib++) { |
| 115 | + if (ohits[ia][ib] >= 1.0) { |
| 116 | + oimage[ia][ib] /= ohits[ia][ib]; |
| 117 | + osqimg[ia][ib] /= ohits[ia][ib]; |
| 118 | + } |
| 119 | + } |
| 120 | + } |
| 121 | + } |
| 122 | + |
| 123 | + if (cram_gather->outaz) { |
| 124 | + for (ia = 0; ia < cram_gather->na/4; ia++) { /* [0; PI) */ |
| 125 | + for (ib = 0; ib < cram_gather->noa; ib++) { |
| 126 | + /* Image */ |
| 127 | + if (oafile) |
| 128 | + sf_floatwrite (&oimage[ia][cram_gather->nb - ib], 1, oafile); |
| 129 | + /* Illumination map */ |
| 130 | + if (ilfile) |
| 131 | + sf_floatwrite (&ohits[ia][cram_gather->nb - ib], 1, ilfile); |
| 132 | + /* Energy image */ |
| 133 | + if (sqfile) |
| 134 | + sf_floatwrite (&osqimg[ia][cram_gather->nb - ib], 1, sqfile); |
| 135 | + } |
| 136 | + } |
| 137 | + for (ia = 0; ia < cram_gather->na/4; ia++) { /* [PI; 2PI) */ |
| 138 | + /* Image */ |
| 139 | + if (oafile) |
| 140 | + sf_floatwrite (&oimage[ia][cram_gather->nb], cram_gather->noa, oafile); |
| 141 | + /* Illumination map */ |
| 142 | + if (ilfile) |
| 143 | + sf_floatwrite (&ohits[ia][cram_gather->nb], cram_gather->noa, ilfile); |
| 144 | + /* Energy image */ |
| 145 | + if (sqfile) |
| 146 | + sf_floatwrite (&osqimg[ia][cram_gather->nb], cram_gather->noa, sqfile); |
| 147 | + } |
| 148 | + } else { |
| 149 | + memset (cram_gather->stack, 0, sizeof(float)*4*cram_gather->nb); |
| 150 | + memset (cram_gather->energy, 0, sizeof(float)*4*cram_gather->nb); |
| 151 | + memset (cram_gather->illum, 0, sizeof(float)*4*cram_gather->nb); |
| 152 | + for (ia = 0; ia < cram_gather->na/4; ia++) { |
| 153 | + for (ib = 0; ib < cram_gather->noa; ib++) { |
| 154 | + cram_gather->stack[ib] += oimage[ia][cram_gather->nb - ib]; |
| 155 | + cram_gather->energy[ib] += osqimg[ia][cram_gather->nb - ib]; |
| 156 | + cram_gather->illum[ib] += ohits[ia][cram_gather->nb - ib]; |
| 157 | + } |
| 158 | + } |
| 159 | + for (ia = 0; ia < cram_gather->na/4; ia++) { |
| 160 | + for (ib = 1; ib < cram_gather->noa; ib++) { |
| 161 | + cram_gather->stack[ib] += oimage[ia][cram_gather->nb + ib]; |
| 162 | + cram_gather->energy[ib] += osqimg[ia][cram_gather->nb + ib]; |
| 163 | + cram_gather->illum[ib] += ohits[ia][cram_gather->nb + ib]; |
| 164 | + } |
| 165 | + } |
| 166 | + /* Image */ |
| 167 | + if (oafile) |
| 168 | + sf_floatwrite (cram_gather->stack, cram_gather->noa, oafile); |
| 169 | + /* Illumination map */ |
| 170 | + if (ilfile) |
| 171 | + sf_floatwrite (cram_gather->illum, cram_gather->noa, ilfile); |
| 172 | + /* Energy image */ |
| 173 | + if (sqfile) |
| 174 | + sf_floatwrite (cram_gather->energy, cram_gather->noa, sqfile); |
| 175 | + } |
| 176 | +} |
| 177 | + |
| 178 | +void sf_cram_gather3_dangle_output (sf_cram_gather3 cram_gather, sf_cram_point3 cram_point, |
| 179 | + sf_file dafile, sf_file sqfile, sf_file ilfile) |
| 180 | +/*< Add one spatial point to the dip angle gather and write it to disk >*/ |
| 181 | +{ |
| 182 | + int ib, ia; |
| 183 | + float **dimage, **dsqimg, **dhits; |
| 184 | + |
| 185 | + dimage = sf_cram_point3_get_dimage (cram_point, &dsqimg, &dhits); |
| 186 | + |
| 187 | + if (cram_gather->inorm) { |
| 188 | + for (ia = 0; ia < cram_gather->na/4; ia++) { |
| 189 | + for (ib = 0; ib < 4*cram_gather->nb; ib++) { |
| 190 | + if (dhits[ia][ib] >= 1.0) { |
| 191 | + dimage[ia][ib] /= dhits[ia][ib]; |
| 192 | + dsqimg[ia][ib] /= dhits[ia][ib]; |
| 193 | + } |
| 194 | + } |
| 195 | + } |
| 196 | + } |
| 197 | + |
| 198 | + if (cram_gather->outaz) { |
| 199 | + for (ia = 0; ia < cram_gather->na/4; ia++) { |
| 200 | + /* Image */ |
| 201 | + if (dafile) |
| 202 | + sf_floatwrite (&dimage[ia][cram_gather->oda], cram_gather->nda, dafile); |
| 203 | + /* Illumination map */ |
| 204 | + if (ilfile) |
| 205 | + sf_floatwrite (&dhits[ia][cram_gather->oda], cram_gather->nda, ilfile); |
| 206 | + /* Energy image */ |
| 207 | + if (sqfile) |
| 208 | + sf_floatwrite (&dsqimg[ia][cram_gather->oda], cram_gather->nda, sqfile); |
| 209 | + } |
| 210 | + } else { |
| 211 | + memset (cram_gather->stack, 0, sizeof(float)*4*cram_gather->nb); |
| 212 | + memset (cram_gather->energy, 0, sizeof(float)*4*cram_gather->nb); |
| 213 | + memset (cram_gather->illum, 0, sizeof(float)*4*cram_gather->nb); |
| 214 | + for (ia = 0; ia < cram_gather->na/4; ia++) { |
| 215 | + for (ib = 0; ib < cram_gather->nda; ib++) { |
| 216 | + cram_gather->stack[ib] += dimage[ia][cram_gather->oda + ib]; |
| 217 | + cram_gather->energy[ib] += dsqimg[ia][cram_gather->oda + ib]; |
| 218 | + cram_gather->illum[ib] += dhits[ia][cram_gather->oda + ib]; |
| 219 | + } |
| 220 | + } |
| 221 | + /* Image */ |
| 222 | + if (dafile) |
| 223 | + sf_floatwrite (cram_gather->stack, cram_gather->nda, dafile); |
| 224 | + /* Illumination map */ |
| 225 | + if (ilfile) |
| 226 | + sf_floatwrite (cram_gather->illum, cram_gather->nda, ilfile); |
| 227 | + /* Energy image */ |
| 228 | + if (sqfile) |
| 229 | + sf_floatwrite (cram_gather->energy, cram_gather->nda, sqfile); |
| 230 | + } |
| 231 | +} |
| 232 | + |
| 233 | +void sf_cram_gather3_setup_image_output (sf_cram_gather3 cram_gather, sf_file input, |
| 234 | + sf_file imfile) |
| 235 | +/*< Setup dimensions in a file for the image >*/ |
| 236 | +{ |
| 237 | + sf_oaxa (imfile, sf_iaxa (input, 4), 1); |
| 238 | + sf_oaxa (imfile, sf_iaxa (input, 5), 2); |
| 239 | + sf_oaxa (imfile, sf_iaxa (input, 6), 3); |
| 240 | + sf_putint (imfile, "n4", 1); |
| 241 | + sf_putint (imfile, "n5", 1); |
| 242 | + sf_putint (imfile, "n6", 1); |
| 243 | +} |
| 244 | + |
| 245 | +void sf_cram_gather3_setup_oangle_output (sf_cram_gather3 cram_gather, sf_file input, |
| 246 | + sf_file oafile) |
| 247 | +/*< Setup dimensions in a file for scattering angle gathers >*/ |
| 248 | +{ |
| 249 | + if (cram_gather->outaz) |
| 250 | + sf_unshiftdim (input, oafile, 0); |
| 251 | + else |
| 252 | + sf_unshiftdim2 (input, oafile, 0); |
| 253 | + sf_putint (oafile, "n1", cram_gather->noa); |
| 254 | + sf_putfloat (oafile, "d1", cram_gather->db); |
| 255 | + sf_putfloat (oafile, "o1", 0.0); |
| 256 | + sf_putstring (oafile, "label1", "Scattering angle"); |
| 257 | + if (cram_gather->outaz) { |
| 258 | + sf_putint (oafile, "n2", cram_gather->na/2); |
| 259 | + sf_putfloat (oafile, "d2", 2.0*cram_gather->da); |
| 260 | + sf_putfloat (oafile, "o2", 0.0); |
| 261 | + sf_putstring (oafile, "label2", "Scattering angle azimuth"); |
| 262 | + } |
| 263 | +} |
| 264 | + |
| 265 | +void sf_cram_gather3_setup_dangle_output (sf_cram_gather3 cram_gather, sf_file input, |
| 266 | + sf_file dafile) |
| 267 | +/*< Setup dimensions in a file for dip angle gathers >*/ |
| 268 | +{ |
| 269 | + if (cram_gather->outaz) |
| 270 | + sf_unshiftdim (input, dafile, 0); |
| 271 | + else |
| 272 | + sf_unshiftdim2 (input, dafile, 0); |
| 273 | + sf_putint (dafile, "n1", cram_gather->nda); |
| 274 | + sf_putfloat (dafile, "d1", 0.5*cram_gather->db); |
| 275 | + sf_putfloat (dafile, "o1", -180.0 + cram_gather->oda*0.5*cram_gather->db); |
| 276 | + sf_putstring (dafile, "label1", "Dip angle"); |
| 277 | + if (cram_gather->outaz) { |
| 278 | + sf_putint (dafile, "n2", cram_gather->na/4); |
| 279 | + sf_putfloat (dafile, "d2", 2.0*cram_gather->da); |
| 280 | + sf_putfloat (dafile, "o2", 0.0); |
| 281 | + sf_putstring (dafile, "label2", "Dip angle azimuth"); |
| 282 | + } |
| 283 | +} |
| 284 | + |
0 commit comments