psvreflect.c

su 的源代码库

							ux[ip]=cadd(cmul(rtb211,vuzs1),cmul(rtb212,vuzs2));
							uz[ip]=cadd(cmul(rtb221,vuzs1),cmul(rtb222,vuzs2));
						}
					}
				} else {
		
					/* receiver below the source */
					if (flag[iz]==2) {
						if (acoustic[iz]==3) {
							for (ip=0; ip<ibl[iz]-1; ip++) {
								ijk1=ik1+ip;
								rtb111=csub(rnr11[ip],rdfr11[ijk1]);
								rtb112=csub(rnr12[ip],rdfr12[ijk1]);
								rtb121=csub(rnr21[ip],rdfr21[ijk1]);
								rtb122=csub(rnr22[ip],rdfr22[ijk1]);
								det=csub(cmul(rtb111,rtb122),
									cmul(rtb112,rtb121));
								rtb211=cdiv(rtb122,det);
								rtb212=cneg(cdiv(rtb112,det));
								rtb221=cneg(cdiv(rtb121,det));
								rtb222=cdiv(rtb111,det);
								rtb311=cadd(cmul(turf11[ijk1],rtb211),
									cmul(rtb212,turf21[ijk1]));
								rtb312=cadd(cmul(turf12[ijk1],rtb211),
									cmul(rtb212,turf22[ijk1]));
								rtb321=cadd(cmul(turf11[ijk1],rtb221),
									cmul(rtb222,turf21[ijk1]));
								rtb322=cadd(cmul(turf12[ijk1],rtb221),
									cmul(rtb222,turf22[ijk1]));
								rtb111=cadd(cmul(tdfr11[ijk1],rtb311),cadd(
									cmul(rtb321,tdfr12[ijk1]),rurf11[ijk1])); 
								rtb112=cadd(cmul(tdfr11[ijk1],rtb312),cadd(
									cmul(rtb322,tdfr12[ijk1]),rurf12[ijk1]));
								rtb121=cadd(cmul(tdfr21[ijk1],rtb311),cadd(
									cmul(rtb321,tdfr22[ijk1]),rurf21[ijk1]));
								rtb122=cadd(cmul(tdfr21[ijk1],rtb312),cadd(
									cmul(rtb322,tdfr22[ijk1]),rurf22[ijk1]));
								det=csub(cmul(rtb111,rtb122),
									cmul(rtb112,rtb121));
								rdnr11[ip]=cdiv(rtb122,det);
								rdnr12[ip]=cneg(cdiv(rtb112,det));
								rdnr21[ip]=cneg(cdiv(rtb121,det));
								rdnr22[ip]=cdiv(rtb111,det);
								rtb111=cadd(cmul(rusr11[ip],rdfr11[ijk1]),
									cmul(rusr12[ip],rdfr21[ijk1]));
								rtb112=cadd(cmul(rusr11[ip],rdfr12[ijk1]),
									cmul(rusr12[ip],rdfr22[ijk1]));
								rtb121=cadd(cmul(rusr21[ip],rdfr11[ijk1]),
									cmul(rusr22[ip],rdfr21[ijk1]));
								rtb122=cadd(cmul(rusr21[ip],rdfr12[ijk1]),
									cmul(rusr22[ip],rdfr22[ijk1]));
								x=csub(cmplx(1.0,0.0),rtb111);
								y=csub(cmplx(1.0,0.0),rtb122);
								det=csub(cmul(x,y),cmul(rtb112,rtb121));
								rvrb111=cdiv(y,det);
								rvrb112=cdiv(rtb112,det);
								rvrb121=cdiv(rtb121,det);
								rvrb122=cdiv(x,det);
								rtb211=cadd(cmul(tdfr11[ijk1],rvrb111),
									cmul(tdfr12[ijk1],rvrb121));
								rtb212=cadd(cmul(tdfr11[ijk1],rvrb112),
									cmul(tdfr12[ijk1],rvrb122));
								rtb221=cadd(cmul(tdfr21[ijk1],rvrb111),
									cmul(tdfr22[ijk1],rvrb121));
								rtb222=cadd(cmul(tdfr21[ijk1],rvrb112),
									cmul(tdfr22[ijk1],rvrb122));
								tdfr11[ijk1]=cadd(cmul(rtb211,tdsr11[ip]),
									cmul(rtb212,tdsr21[ip]));
								tdfr12[ijk1]=cadd(cmul(rtb211,tdsr12[ip]),
									cmul(rtb212,tdsr22[ip]));
								tdfr21[ijk1]=cadd(cmul(rtb221,tdsr11[ip]),
									cmul(rtb222,tdsr21[ip]));
								tdfr22[ijk1]=cadd(cmul(rtb221,tdsr12[ip]),
									cmul(rtb222,tdsr22[ip]));
								rtb111=cadd(cmul(rusr11[ip],rtb211),cmul(rtb212,
									rusr21[ip]));
								rtb112=cadd(cmul(rusr12[ip],rtb211),cmul(rtb212,
									rusr22[ip]));
								rtb121=cadd(cmul(rusr11[ip],rtb221),cmul(rtb222,
									rusr21[ip]));
								rtb122=cadd(cmul(rusr12[ip],rtb221),cmul(rtb222,
									rusr22[ip]));
								rurf11[ijk1]=cadd(cmul(rtb111,rurf11[ijk1]),
									cadd(cmul(rtb112,turf21[ijk1]),
									rurf11[ijk1]));
								rurf12[ijk1]=cadd(cmul(rtb111,rurf12[ijk1]),
									cadd(cmul(rtb112,turf22[ijk1]),
									rurf12[ijk1]));
								 rurf21[ijk1]=cadd(cmul(rtb121,rurf11[ijk1]),
									cadd(cmul(rtb122,turf21[ijk1]),
									rurf21[ijk1]));
								rurf22[ijk1]=cadd(cmul(rtb121,rurf12[ijk1]),
									cadd(cmul(rtb122,turf22[ijk1]),
									rurf22[ijk1]));
							}
							for (ip=ibl[iz]-1; ip<block_size; ip++) {
								ijk1=ik1+ip;
								rdnr11[ip]=cmplx(0.0,0.0);
								rdnr12[ip]=cmplx(0.0,0.0);
								rdnr21[ip]=cmplx(0.0,0.0);
								rdnr22[ip]=cmplx(0.0,0.0);
								rtb111=cadd(cmul(rusr11[ip],rdfr11[ijk1]),
									cmul(rusr12[ip],rdfr21[ijk1]));
								rtb112=cadd(cmul(rusr11[ip],rdfr12[ijk1]),
									cmul(rusr12[ip],rdfr22[ijk1]));
								rtb121=cadd(cmul(rusr21[ip],rdfr11[ijk1]),
									cmul(rusr22[ip],rdfr21[ijk1]));
								rtb122=cadd(cmul(rusr21[ip],rdfr12[ijk1]),
									cmul(rusr22[ip],rdfr22[ijk1]));
								x=csub(cmplx(1.0,0.0),rtb111);
								y=csub(cmplx(1.0,0.0),rtb122);
								det=csub(cmul(x,y),cmul(rtb112,rtb121));
								rvrb111=cdiv(y,det);
								rvrb112=cdiv(rtb112,det);
								rvrb121=cdiv(rtb121,det);
								rvrb122=cdiv(x,det);
								rtb211=cadd(cmul(tdfr11[ijk1],rvrb111),
									cmul(rvrb121,tdfr12[ijk1]));
								rtb212=cadd(cmul(tdfr11[ijk1],rvrb112),
									cmul(rvrb122,tdfr12[ijk1]));
								rtb221=cadd(cmul(tdfr21[ijk1],rvrb111),
									cmul(rvrb121,tdfr22[ijk1]));
								rtb222=cadd(cmul(tdfr21[ijk1],rvrb112),
									cmul(rvrb122,tdfr22[ijk1]));
								tdfr11[ijk1]=cadd(cmul(rtb211,tdsr11[ip]),	
									cmul(rtb212,tdsr21[ip]));
								tdfr12[ijk1]=cadd(cmul(rtb211,tdsr12[ip]),
									cmul(rtb212,tdsr22[ip]));
								tdfr21[ijk1]=cadd(cmul(rtb221,tdsr11[ip]),
									cmul(rtb222,tdsr21[ip]));
								tdfr22[ijk1]=cadd(cmul(rtb221,tdsr12[ip]),
									cmul(rtb222,tdsr22[ip]));
								rtb111=cadd(cmul(rtb211,rusr11[ip]),
									cmul(rtb212,rusr21[ip]));
								rtb112=cadd(cmul(rtb211,rusr12[ip]),
									cmul(rtb212,rusr22[ip]));
								rtb121=cadd(cmul(rtb221,rusr11[ip]),
									cmul(rtb222,rusr21[ip]));
								rtb122=cadd(cmul(rtb221,rusr12[ip]),
									cmul(rtb222,rusr22[ip]));
								rurf11[ijk1]=cadd(cmul(rtb111,turf11[ijk1]),
									cadd(cmul(rtb112,turf21[ijk1]),
									rurf11[ijk1]));
								rurf12[ijk1]=cadd(cmul(rtb111,turf12[ijk1]),
									cadd(cmul(rtb112,turf22[ijk1]),
									rurf12[ijk1]));
								rurf21[ijk1]=cadd(cmul(rtb121,turf11[ijk1]),
									cadd(cmul(rtb122,turf21[ijk1]),
								rurf21[ijk1]));
								rurf22[ijk1]=cadd(cmul(rtb121,turf12[ijk1]),
									cadd(cmul(rtb122,turf22[ijk1]),
									rurf22[ijk1]));
							}
						} else if (acoustic[iz]==1) {
							for (ip=0; ip<ibl[iz]-1; ip++) {
								ijk1=ik1+ip;
								rtb211=cdiv(turf11[ijk1],csub(rd0n11[ip],
									rdfr11[ijk1]));
								rdnr11[ip]=cdiv(cmplx(1.0,0.0),cadd(cmul(tdfr11[ijk1],
									rtb211),rurf11[ijk1]));
								rdnr12[ip]=cmplx(0.0,0.0);   
								rdnr21[ip]=cmplx(0.0,0.0);   
								rdnr22[ip]=cmplx(0.0,0.0);   
							}
							for (ip=ibl[iz]-1; ip<block_size; ip++) {
								rdnr11[ip]=cmplx(0.0,0.0);
								rdnr12[ip]=cmplx(0.0,0.0);
								rdnr21[ip]=cmplx(0.0,0.0);
								rdnr22[ip]=cmplx(0.0,0.0);
							}
						} else if (acoustic[iz]==2) {
							for (ip=0; ip<block_size; ip++) {
								rdnr11[ip]=rnr11[ip];
								rdnr12[ip]=rnr12[ip];
								rdnr21[ip]=rnr21[ip];
								rdnr22[ip]=rnr22[ip];
							}
						} else {
							for (ip=0; ip<ibl[iz]-1; ip++) {
								ijk1=ik1+ip;
								rtb111=csub(rd0n11[ip],rdfr11[ijk1]);
								rtb112=csub(rd0n12[ip],rdfr12[ijk1]);
								rtb121=csub(rd0n21[ip],rdfr21[ijk1]);
								rtb122=csub(rd0n22[ip],rdfr22[ijk1]);
								det=csub(cmul(rtb111,rtb122),
									cmul(rtb112,rtb121));
								rtb211=cdiv(rtb122,det);
								rtb212=cneg(cdiv(rtb112,det));
								rtb221=cneg(cdiv(rtb121,det));
								rtb222=cdiv(rtb111,det);
								rtb311=cadd(cmul(rtb211,turf11[ijk1]),
									cmul(rtb212,turf21[ijk1]));
								rtb312=cadd(cmul(rtb211,turf12[ijk1]),
									cmul(rtb212,turf22[ijk1]));
								rtb321=cadd(cmul(rtb221,turf11[ijk1]),
									cmul(rtb222,turf21[ijk1]));
								rtb322=cadd(cmul(rtb221,turf12[ijk1]),
									cmul(rtb222,turf22[ijk1]));
								rtb111=cadd(cmul(rtb311,tdfr11[ijk1]),cadd(
									cmul(rtb321,tdfr12[ijk1]),rurf11[ijk1]));
								rtb112=cadd(cmul(rtb312,tdfr11[ijk1]),cadd(
									cmul(rtb322,tdfr12[ijk1]),rurf12[ijk1]));
								rtb121=cadd(cmul(rtb311,tdfr21[ijk1]),cadd(
									cmul(rtb321,tdfr22[ijk1]),rurf21[ijk1]));
								rtb122=cadd(cmul(rtb312,tdfr21[ijk1]),cadd(
									cmul(rtb322,tdfr22[ijk1]),rurf22[ijk1]));
								det=csub(cmul(rtb111,rtb122),
									cmul(rtb112,rtb121));
								rdnr11[ip]=cdiv(rtb122,det);
								rdnr12[ip]=cneg(cdiv(rtb112,det));
								rdnr21[ip]=cneg(cdiv(rtb121,det));
								rdnr22[ip]=cdiv(rtb111,det);
							}
							for (ip=ibl[iz]-1; ip<block_size; ip++) {
								rdnr11[ip]=cmplx(0.0,0.0);
								rdnr12[ip]=cmplx(0.0,0.0);
								rdnr21[ip]=cmplx(0.0,0.0);
								rdnr22[ip]=cmplx(0.0,0.0);
							}
						}
					} else if (flag[iz]==0) {
						/* source acoustic, 1st elastic rec */
						if (acoustic[iz]==2) {
							for (ip=0; ip<block_size; ip++) {
								rdnr11[ip]=rnr11[ip];
								rdnr12[ip]=rnr12[ip];
								rdnr21[ip]=rnr21[ip];
								rdnr22[ip]=rnr22[ip];
							}
						} else if (acoustic[iz]==3) {
							for (ip=0; ip<block_size; ip++) {
								ijk1=ik1+ip;
								rtb111=csub(rnr11[ip],rdfr11[ijk1]);
								rtb112=csub(rnr12[ip],rdfr12[ijk1]);
								rtb121=csub(rnr21[ip],rdfr21[ijk1]);
								rtb122=csub(rnr22[ip],rdfr22[ijk1]);
								det=csub(cmul(rtb111,rtb122),
									cmul(rtb112,rtb121));
								rtb211=cdiv(rtb122,det);
								rtb212=cneg(cdiv(rtb112,det));
								rtb221=cneg(cdiv(rtb121,det));
								rtb222=cdiv(rtb111,det);
								rtb311=cadd(cmul(turf11[ijk1],rtb211),
									cmul(rtb212,turf21[ijk1]));
								rtb312=cadd(cmul(turf12[ijk1],rtb211),
									cmul(rtb212,turf22[ijk1]));
								rtb321=cadd(cmul(turf11[ijk1],rtb221),
									cmul(rtb222,turf21[ijk1]));
								rtb322=cadd(cmul(turf12[ijk1],rtb221),
									cmul(rtb222,turf22[ijk1]));
								rtb111=cadd(cmul(tdfr11[ijk1],rtb311),cadd(
									cmul(rtb321,tdfr12[ijk1]),rurf11[ijk1]));
								rtb112=cadd(cmul(tdfr11[ijk1],rtb312),cadd(
									cmul(rtb322,tdfr12[ijk1]),rurf12[ijk1]));
								rtb121=cadd(cmul(tdfr21[ijk1],rtb311),cadd(
									cmul(rtb321,tdfr22[ijk1]),rurf21[ijk1]));
								rtb122=cadd(cmul(tdfr21[ijk1],rtb312),cadd(
									cmul(rtb322,tdfr22[ijk1]),rurf22[ijk1]));
								det=csub(cmul(rtb111,rtb122),
									cmul(rtb112,rtb121));
								rdnr11[ip]=cdiv(rtb122,det);
								rdnr12[ip]=cneg(cdiv(rtb112,det));
								rdnr21[ip]=cneg(cdiv(rtb121,det));
								rdnr22[ip]=cdiv(rtb111,det);
								rtb111=cadd(cmul(rusr11[ip],rdfr11[ijk1]),
									cmul(rusr12[ip],rdfr21[ijk1]));
								rtb112=cadd(cmul(rusr11[ip],rdfr12[ijk1]),
									cmul(rusr12[ip],rdfr22[ijk1]));
								rtb121=cadd(cmul(rusr21[ip],rdfr11[ijk1]),
									cmul(rusr22[ip],rdfr21[ijk1]));
								rtb122=cadd(cmul(rusr21[ip],rdfr12[ijk1]),
									cmul(rusr22[ip],rdfr22[ijk1]));
								x=csub(cmplx(1.0,0.0),rtb111);
								y=csub(cmplx(1.0,0.0),rtb122);
								det=csub(cmul(x,y),cmul(rtb112,rtb121));
								rvrb111=cdiv(y,det);
								rvrb112=cdiv(rtb112,det);
								rvrb121=cdiv(rtb121,det);
								rvrb122=cdiv(x,det);
								rtb211=cadd(cmul(tdfr11[ijk1],rvrb111),
									cmul(rvrb121,tdfr12[ijk1]));
								rtb212=cadd(cmul(tdfr11[ijk1],rvrb112),
									cmul(rvrb122,tdfr12[ijk1]));
								rtb221=cadd(cmul(tdfr21[ijk1],rvrb111),
									cmul(rvrb121,tdfr22[ijk1]));
								rtb222=cadd(cmul(tdfr21[ijk1],rvrb112),
									cmul(rvrb122,tdfr22[ijk1]));
								tdfr11[ijk1]=cadd(cmul(rtb211,tdsr11[ip]),
									cmul(rtb212,tdsr21[ip]));
								tdfr12[ijk1]=cadd(cmul(rtb211,tdsr12[ip]),
									cmul(rtb212,tdsr22[ip]));
								tdfr21[ijk1]=cadd(cmul(rtb221,tdsr11[ip]),
									cmul(rtb222,tdsr21[ip]));
								tdfr22[ijk1]=cadd(cmul(rtb221,tdsr12[ip]),
									cmul(rtb222,tdsr22[ip]));
								rtb111=cadd(cmul(rtb211,rusr11[ip]),
									cmul(rtb212,rusr21[ip]));
								rtb112=cadd(cmul(rtb211,rusr12[ip]),cmul(rtb212,
									rusr22[ip]));
								rtb121=cadd(cmul(rtb221,rusr11[ip]),cmul(rtb222,
									rusr21[ip]));
								rtb122=cadd(cmul(rtb221,rusr12[ip]),cmul(rtb222,
									rusr22[ip]));
								rurf11[ijk1]=cadd(cmul(rtb111,turf11[ijk1]),
									cadd(cmul(rtb112,turf21[ijk1]),
									rurf11[ijk1]));
								rurf12[ijk1]=cadd(cmul(rtb111,turf12[ijk1]),
									cadd(cmul(rtb112,turf22[ijk1]),
									rurf12[ijk1]));
								rurf21[ijk1]=cadd(cmul(rtb121,turf11[ijk1]),
									cadd(cmul(rtb122,turf21[ijk1]),
									rurf21[ijk1]));
								rurf22[ijk1]=cadd(cmul(rtb121,turf12[ijk1]),
									cadd(cmul(rtb122,turf22[ijk1]),
									rurf22[ijk1]));
							}
						} else if (acoustic[iz]==0) {
							for (ip=0; ip<block_size; ip++) {
								ijk1=ik1+ip;
								rtb111=csub(rd0n11[ip],rdfr11[ijk1]);
								rtb112=csub(rd0n12[ip],rdfr12[ijk1]);
								rtb121=csub(rd0n21[ip],rdfr21[ijk1]);
								rtb122=csub(rd0n22[ip],rdfr22[ijk1]);
								det=csub(cmul(rtb111,rtb122),
									cmul(rtb112,rtb121));
								rtb211=cdiv(rtb122,det);
								rtb212=cneg(cdiv(rtb112,det));
								rtb221=cneg(cdiv(rtb121,det));
								rtb222=cdiv(rtb111,det);
								rtb311=cadd(cmul(rtb211,turf11[ijk1]),
									cmul(rtb212,turf21[ijk1]));
								rtb312=cadd(cmul(rtb211,turf12[ijk1]),
									cmul(rtb212,turf22[ijk1]));
								rtb321=cadd(cmul(rtb221,turf11[ijk1]),
									cmul(rtb222,turf21[ijk1]));
								rtb322=cadd(cmul(rtb221,turf12[ijk1]),	
									cmul(rtb222,turf22[ijk1]));
								rtb111=cadd(cmul(tdfr11[ijk1],rtb311),cadd(
									cmul(tdfr12[ijk1],rtb321),rurf11[ijk1]));
								rtb112=cadd(cmul(tdfr11[ijk1],rtb312),cadd(
									cmul(tdfr12[ijk1],rtb322),rurf12[ijk1]));
								rtb121=cadd(cmul(tdfr21[ijk1],rtb311),cadd(
									cmul(tdfr22[ijk1],rtb321),rurf21[ijk1]));
								rtb122=cadd(cmul(tdfr21[ijk1],rtb312),cadd(
									cmul(tdfr22[ijk1],rtb322),rurf22[ijk1]));
								det=csub(cmul(rtb111,rtb122),
									cmul(rtb112,rtb121));
								rdnr11[ip]=cdiv(rtb122,det);
								rdnr12[ip]=cneg(cdiv(rtb112,det));
								rdnr21[ip]=cneg(cdiv(rtb121,det));
								rdnr22[ip]=cdiv(rtb111,det);
							}
						} else if (acoustic[iz]==1) {
							for (