PDA

View Full Version : فاصله ی نقاط در شبکه سه بعدی



hamidkor
چهارشنبه 31 اردیبهشت 1393, 19:17 عصر
باسلام دوستان عزیز، کد زیر رو واسه تعیین حالت های مختلف فواصل نقاط نوشته شده است (هرچند ماتریس دادهها بصورت تصادفی تعیین شده است که جهت ارائه مثال می باشد)ولی میخواهیم فاصله نقاطی که نسبت به هم محاسبه شده اند بار دوم محاسبه نشود. به عبارت دیگر اگر فاصله دو نقطه b از a باشد h آنگاه فاصله a تا b نیز h خواهد بود که ما نمیخواهیم هر دو ثبت گردد و تنها یک مورد کافی است. همچنین فاصله ها حداکثر تا نصف محدوده مورد بررسی قابل قبول می باشد. درضمن هرگونه پیشنهاد دوستان برای بهبود سرعت اجرا و سایر موارد ارزشمند خواهد بود.


a=rand(10000,103);
k=0;
xmax=max(a(:,1));
ymax=max(a(:,2));
zmax=max(a(:,3));
xmin=min(a(:,1));
ymin=min(a(:,2));
zmin=min(a(:,3));
hmax=(((xmax-xmin)^2)+((ymax-ymin)^2)+((zmax-zmin)^2))^0.5;
for i=1:size(a,1)-1
for j=i+1:size(a,1)
x=a(i,1)-a(j,1);
y=a(i,2)-a(j,2);
z=a(i,3)-a(j,3);
while ((x^2)+(y^2)+(z^2))^(1/2)<=(hmax/2);
k=k+1;
h(k,1)=((x^2)+(y^2)+(z^2))^(1/2);
h(k,4:end)=(a(i,4:end)-a(j,4:end)).^2;
end
end
end

rahnema1
چهارشنبه 31 اردیبهشت 1393, 21:52 عصر
سلام
فکر کنم این جور بشه

order = nchoosek(1:rows(a),2);
Xi = order(:,1);
Yi = order(:,2);
d = a(Xi,1:3) - a(Yi,1:3);
y = sqrt(sumsq (d,2));
lowindex=find(y<=max(y));
h=[y(lowindex) , (a(Xi(lowindex),4:end)-a(Yi(lowindex),4:end)).^2];

hamidkor
پنج شنبه 01 خرداد 1393, 20:05 عصر
سلام
فکر کنم این جور بشه

order = nchoosek(1:rows(a),2);
Xi = order(:,1);
Yi = order(:,2);
d = a(Xi,1:3) - a(Yi,1:3);
y = sqrt(sumsq (d,2));
lowindex=find(y<=max(y));
h=[y(lowindex) , (a(Xi(lowindex),4:end)-a(Yi(lowindex),4:end)).^2];


باسلام،
ضمن تشکر بابت کدتون، لازم دونستم خطای تعریف نشده بودن تابع sumsq رو واسه متلب بیان کنم و اینکه شرط ما تا نصف محدوده هستش و الگوریتم برای یک شبکه ی منظم سه بعدی از نقاط هستش.

hamidkor
پنج شنبه 01 خرداد 1393, 20:06 عصر
سلام
فکر کنم این جور بشه

order = nchoosek(1:rows(a),2);
Xi = order(:,1);
Yi = order(:,2);
d = a(Xi,1:3) - a(Yi,1:3);
y = sqrt(sumsq (d,2));
lowindex=find(y<=max(y));
h=[y(lowindex) , (a(Xi(lowindex),4:end)-a(Yi(lowindex),4:end)).^2];


باسلام،
ضمن تشکر بابت کدتون، لازم دونستم خطای تعریف نشده بودن تابع sumsq رو واسه متلب بیان کنم و اینکه شرط ما تا نصف محدوده هستش و الگوریتم برای یک شبکه ی منظم سه بعدی از نقاط می باشد.

rahnema1
پنج شنبه 01 خرداد 1393, 22:05 عصر
ببخشید من توی octave تست کردم که ظاهرا توی متلب نیست . اینجور فکر کنم حل بشه

order = nchoosek(1:rows(a),2);
Xi = order(:,1);
Yi = order(:,2);
d = a(Xi,1:3) - a(Yi,1:3);
y = sqrt(sum (d.^2,2));
lowindex=find(y<=max(y)/2);
h=[y(lowindex) , (a(Xi(lowindex),4:end)-a(Yi(lowindex),4:end)).^2];

rahnema1
جمعه 02 خرداد 1393, 16:32 عصر
یه اصلاح کردم حواسم نبود تقسیم بر 2 کنم

hamidkor
پنج شنبه 05 تیر 1393, 14:06 عصر
دوست عزیز ضمن تشکر مجدد بابت راهنمایی های ارزندتون میخواستم خدمتتان عارض بشم که بنده کد شما رو و کد خودم رو بصورت زیر برای رسیدن به هدفم ادامه دادم ولی موقع تست از آنجائیکه با داده های کم حجم امتحان کرده بودم از امکان به مشکل خوردن سر داده های حجیم غافل بودم و این شد که الان بابت مشکل جدید که زمان بر بودن فرآیند و در نهایت هنگ کردن سیستم برای داده هایی با ابعاد مثلا 150000 در 100 هستش مجددا مزاحم وقت جنابعالی شده ام.



% a: a matrices with large dimension
order = nchoosek(1:size(a,1),2);
Xi = order(:,1);
Yi = order(:,2);
d = a(Xi,1:3) - a(Yi,1:3);
y = sqrt(sum (d.^2,2));
lowindex=find(y<=(max(y)/2));
h=[y(lowindex) , (a(Xi(lowindex),4:end)-a(Yi(lowindex),4:end)).^2];
u=unique(h(:,1));
c=zeros(size(u,1),size(h,2));
for i=1:numel(u)
c(i,:)=[u(i,1) sum(h(h(:,1)==u(i,1),2:end),1)];
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%
countOfu = hist(h(:,1),u);
b = [u countOfu'];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
g=zeros(size(c));
for i=1:size(b,1)
g(i,:)=[c(i,1) (c(i,2:end))./(2*b(i,2))];
end

hamidkor
پنج شنبه 05 تیر 1393, 15:42 عصر
داده های ماتریس اولیه اگر لازم باشند قابل ارائه می باشند.

hamidkor
پنج شنبه 05 تیر 1393, 15:49 عصر
کد اولیه نیز برای فلوچارت اصلی بصورت زیر اصلاح شده است که برای داده های کم حجم جواب خوبی ارائه می دهد. لازم به یادآوری است که محاسبات برای شبکه منظم سه بعدی از داده ها است.

120513



%a:a matrices with large dimension
k=0;
xmax=max(a(:,1));
ymax=max(a(:,2));
zmax=max(a(:,3));
xmin=min(a(:,1));
ymin=min(a(:,2));
zmin=min(a(:,3));
hmax=(((xmax-xmin)^2)+((ymax-ymin)^2)+((zmax-zmin)^2))^0.5;
h=zeros(size(a,1),size(a,2)-2);
for i=1:size(a,1)-1
for j=i+1:size(a,1)
x=a(i,1)-a(j,1);
y=a(i,2)-a(j,2);
z=a(i,3)-a(j,3);
if ((x^2)+(y^2)+(z^2))^(1/2)<=((hmax)/2);
k=k+1;
h(k,1)=((x^2)+(y^2)+(z^2))^(1/2);
h(k,2:end)=((a(i,4:end)-a(j,4:end)).^2);
end
end
end

u=unique(h(:,1));
c=zeros(size(u,1),size(h,2));
for i=1:numel(u)
c(i,:)=[u(i,1) sum(h(h(:,1)==u(i,1),2:end),1)];
end

%%%%%%%%%%%%%%%%%%%%%%%%%
countOfu = hist(h(:,1),u);
b = [u countOfu'];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
g=zeros(size(c));
for i=1:size(b,1)
g(i,:)=[c(i,1) (c(i,2:end))./(2*b(i,2))];
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
plot(g(:,1),g(:,2:end))
set(gca,'XTick',g(:,1))
title('Realizations Variogram');
xlabel('Distance');
ylabel('Variogram Value');

rahnema1
پنج شنبه 05 تیر 1393, 16:58 عصر
مگه قرار نیست سه بعدی باشه پس چرا داده های دیگه را دخالت می دهید؟ یعنی چه لزومی داره؟

h(k,2:end)=((a(i,4:end)-a(j,4:end)).^2);

یک نکته دیگه حالا میشه یه جوری با حلقه چاره ای براش پیدا کرد اما تنها محاسبه lowindex نزدیگ 4 گیگ حافظه لازم داره
فلوچارت کاملا گویا نیست class یعنی چی؟

hamidkor
پنج شنبه 05 تیر 1393, 21:37 عصر
درسته سه بعدیه و داده های دیگه مقادیر نقاطی هستند که مختصاتشان را در سه بعد داریم و ما نقاطی را که در سه بعد به فاصله های یکسان از هم هستند را در یک کلاس قرار می دهیم و مجذور اختلاف مقادیر (عیار نقاط) آنها نیز با هم جمع میشود و سرآخر به دو برابر تعداد جفت نقاط هر کلاس تقسیم می شود. بصورت فرمولی که در تصویر دیده می شود.
120519
در این رابطه xi موقعیت یامختصات هر نقطه است و z(xi) مقادیر آن نقطه است و N(h) تعداد جفت نقاطی است که با فاصله h ازهم جدا شده اند و h های مختلف کلاس های مختلف را تشکیل می دهند.

ما در این الگوریتم مطابق فلوچارت اومدیم ابتدا فواصل تمام جفت نقاط امکان پذیر و مجذور اختلاف مقادیر متناظر جفت نقاط را محاسبه کردیم و در ادامه به دسته بندی آنها پرداختیم و در نهایت مجموع مقادیر هر کلاس را به دوبرابر تعداد جفت نقاط هر کلاس تقسیم کردیم. احتمالا مرحله ی یافتن فواصل تمام جفت نقاط ممکن سخت ترین و زمان بر ترین مرحله است که به قول شما حافظه بالایی هم نیاز دارد که سیستمی که من در اختیار دارم تنها 3 گیگ رم داره و احتمالا در این مرحله بامشکل مواجه می شه.

rahnema1
پنج شنبه 05 تیر 1393, 21:45 عصر
این جوری نمیشه بهتره از داده ها به صورت رندوم sample بگیرید بعد محاسبه روی اونها انجام بشه. چون ماتریس 1000000 در 1000000 باید تشکیل بشه که عملی نیست

rahnema1
پنج شنبه 05 تیر 1393, 22:36 عصر
یک نکته دیگه فکر کنم فاصله ها اشتباه حساب شده. فکر کنم باید فاصله دوبعدی حساب بشه

hamidkor
جمعه 06 تیر 1393, 07:44 صبح
این جوری نمیشه بهتره از داده ها به صورت رندوم sample بگیرید بعد محاسبه روی اونها انجام بشه. چون ماتریس 1000000 در 1000000 باید تشکیل بشه که عملی نیست

حق با شماست اگر نخواهیم هیچ محدودیتی در نظر بگیریم جفت نقاط مختلف زیادی وجود خواهند داشت. اما اگر این محاسبات برای نصف محدوده انجام بگیرد معتبر خواهد بود و فاصله ی مینیمم را میتوان بزرگتر در نظر گرفت همانطور که در نرم افزارهای صنعتی که برای این کار در نظر گرفته شده اند h را خود کاربر می تواند تعریف کند و کلاس ها مضرب هایی از h تعریف شده خواهند بود. از این دو مورد نصف محدوده در الگوریتم شما در قسمت y<=(max(y)/2)اعمال شده است و در الگوریتم ابتدایی هم در قسمت if و مورد دوم را واقعیتش هنوز بهش فک نکردم که چه جوری می توان اعمال نمود. اما در مورد نصف محدوده اگر بتوان در الگوریتم شما در مراحل اولیه اعمال نمود خیلی از حجم محاسبات کم خواهد شد. همچنین مورد مینیمم h و مضارب آن نیز می توان خیلی تو کاهش حجم محاسبات موثر باشد البته اگر بتوان اعمال نمود.

hamidkor
جمعه 06 تیر 1393, 07:49 صبح
یک نکته دیگه فکر کنم فاصله ها اشتباه حساب شده. فکر کنم باید فاصله دوبعدی حساب بشه

درسته دوست عزیز با توجه به فلوچارت و فرمول ارائه شده که هر دو برای فواصل در دو بعد هستند گفته شما کاملا صحیحه. ولی فلوچارت و فرمول برای رساندن مفهوم کلی الگوریتم ارائه شده و الگوریتم در واقع تعمیم یافته فلوچارته است و از آنجائیکه هیچ فلوچارتی برای حالت سه بعدی یافت نشد و همچنین فرمول نیز می تواند برای حالت سه بعدی تعمیم یابد (چون مهم یافتن دو نقطه ای است که فاصله ی hاز هم داشته باشند و مطابق فرمول فوق مقادیر واریوگرام برای آنها محاسبه گردد) بنابراین به ناچار این فلوچارت و فرمول خدمتتان ارائه گردید.

rahnema1
جمعه 06 تیر 1393, 22:15 عصر
برای سرعت بیشتر بهتره در ++c نوشته می شد که فایلش را قرار دادم
فرض کردم که فایل داده شما مثلا به نام a.txt که شامل 1000000 سطر و 103 ستون هست و داده ها با space از هم جدا شدند
این فایل را می تونید با dlmwrite به این صورت تولید کنید

dlmwrite('a.txt',a,' ')

برای اجرا برنامه لازمه که فایل a.txt در همون فولدری باشه که فایل exe قرار داره همچنین در برنامه تعداد کلاسهای یا bin ها را درخواست می کنه که باید وارد کنید نهایتا خروجی به صورت یک فایل به نام b.txt هست که ستون اول کلاسهای فاصله و ستونهای بعدی هم که 100 ستون بعدی هست و تعداد سطرها هم برابر تعداد bin هست
http://www.sharefile.ir/uploads/1403954319.zip

hamidkor
جمعه 06 تیر 1393, 23:19 عصر
برای سرعت بیشتر بهتره در ++c نوشته می شد که فایلش را قرار دادم
فرض کردم که فایل داده شما مثلا به نام a.txt که شامل 1000000 سطر و 103 ستون هست و داده ها با space از هم جدا شدند
این فایل را می تونید با dlmwrite به این صورت تولید کنید

dlmwrite('a.txt',a,' ')

برای اجرا برنامه لازمه که فایل a.txt در همون فولدری باشه که فایل exe قرار داره همچنین در برنامه تعداد کلاسهای یا bin ها را درخواست می کنه که باید وارد کنید نهایتا خروجی به صورت یک فایل به نام b.txt هست که ستون اول کلاسهای فاصله و ستونهای بعدی هم که 100 ستون بعدی هست و تعداد سطرها هم برابر تعداد bin هست
http://www.sharefile.ir/uploads/1403954319.zip

ممنون دوست عزیز مثل همیشه گل کاشتین، فقط یه سوال: تعداد سطرها اگر کمتر یا بیشتر از این مقدار باشند لازمه که در برنامه تغییری داده بشه ویا اینکه از داده های ورودی برایش قابل تعیین هست؟

rahnema1
شنبه 07 تیر 1393, 03:32 صبح
توی این یکی قابل تنظیم هست
http://www.sharefile.ir/uploads/1404006458.zip

hamidkor
شنبه 07 تیر 1393, 16:23 عصر
توی این یکی قابل تنظیم هست
http://www.sharefile.ir/uploads/1404006458.zip

اطلاعات درخواستی رو وارد میکنم و کلی هم منتظر می مونم (واسه بار اول چند ساعت و واسه کد فایل آخری با وارد کردن سه پارامتر درخواستیش نزدیک یه ساعت) ولی جوابی نمیده و تنها با CPU بالا (100 درصد) سیستم کار میکنه ...

بنظر شما تو کدوم قسمت به مشکل می خوره؟ تو ساخت ماتریس order (order = nchoosek(1:size(a,1),2);) یا تو قسمت find و یا اینکه تو حلقه ی for تودرتو؟

rahnema1
شنبه 07 تیر 1393, 17:00 عصر
این کدی که می فرستم مقدار زمانی که محاسبه برای 1000000 رکورد مصرف خواهد کرد ( به ساعت) را چاپ می کنه لطفا فایل a.txt در فولدر فایل exe باشه
علتش اینه که تعداد محاسبات زیاد هست
http://www.sharefile.ir/uploads/1404060128.zip

rahnema1
یک شنبه 08 تیر 1393, 06:59 صبح
اون برنامه ای که زمان را حساب می کنه یه اشتباه کردم صفرها را زیاد گذاشتم که اصلاح کردم اما باز هم محاسبه 1000000 رکورد زمان خیلی زیادی می گیره
http://www.sharefile.ir/uploads/1404047299.zip
پیشنهاد من اینه که شما با استفاده از randperm از سطرهای ماتریس مثلا به تعداد 10000 تا نمونه بگیرید و واریو گرام را برای این داده های نمونه گیری شده بدست بیارید. حتی می تونید کار نمونه گیری را مثلا 5 بار انجام بدید و 5 تا فایل b.txt درست کنید. فکر می کنم جواب هر فایل شبیه هم باشه. حتما لازم نیست محاسبه روی کل داده ها انجام بشه
در صورت نیاز هم می تونید از اون 5 تا فایل میانگین بگیرید که خاطر جمع بشید

hamidkor
سه شنبه 10 تیر 1393, 18:33 عصر
سلام دوست عزیز

ممنون از این همه لطف و توجهتون، واقعیتش اینکه تمام فاصله ها رو حساب کنیم لازم نیست و من دنبال اینم که نصف محدوده را قبل از محاسبه فاصله ها اعمال کنم و نرم افزارای این کار هم قبلنا میومدن این کار رو رو سه صفحه اصلی مختصات انجام میدادن و نصف فاصله هم قاعدتااونجا اعمال میشده. ولی جدیدا یکی ازنرم افزارای این کار که خیلی هم گرونه و خارجیه همه رو در مدت زمان کمی ارائه میده شاید تو سه صفحه محاسبه و میانگین میگیره احتمالا! اینکه سمپل گرفته بشه هم فکر خیلی جالبیه و اگر واریوگرام حاصل با واریوگرام اولیه شبیه یا نزدیک به هم بشه مشکل حل خواهد شد. فقط یه موضوعی که پیش میاد اینکه بعد سمپل گیری تصادفی آیا باید با کد سی پلاپلاس اجرا بشه و یا میشه با کد اولی که شما تو متلب لطف کردین کار رو راانداخت؟

rahnema1
سه شنبه 10 تیر 1393, 19:19 عصر
من بعید می دونم یک نرم افزار بیاد روی تمام داده ها محاسبه بکنه و مدت زمان کمی هم مصرف بشه. حتی اگه فرض کنیم نرم افزار بیاد محاسبات را روی هسته های سی پی یو هم تقسیم کنه مثلا سی پی یوی شما 4 تا هم هسته داشته باشه فوقش محاسبات یک چهارم خواهد شد. در ضمن ما داریم برای 100 تا ستون محاسبه می کنیم که بیشتر زمان را این مورد مصرف می کنه. اما اگه محاسبه تنها برای یک ستون بود طبیعتا زمان کمتری مصرف می شد.
من فکر می کنم نرم افزار ها سمپیل می گیرند. به هر حال اونها هم با سی پلاس پلاس نوشته میشن و کار جادویی انجام نمی دن
اگه بخواهیم با متلب همین را انجام بدیم فکر کنم زمان بیشتری می گیره و اون کدی که گذاشته بودم واریوگرام را محاسبه نمی کرد. فکر کنم با همون فایل exe بهتر باشه

hamidkor
سه شنبه 10 تیر 1393, 19:44 عصر
راسیتش سمپلشو نمیدونم که میگیرن یا نه ولی از آنجائیکه خیلی خیلی گرونه و نمیشه لایسنسشو تهیه کرد و از طرفی قفل شکسته هاشم خروجی دیتا نمیدن و تنها نمودارشو رو میکنه و تنها میتونه میانگینشون در اختیار بده بنابراین نمیشه نمودارشو مثلا تو متلب یا اکسل هم تشکیل داد و نرم افزاره رو دورش زد بنابراین واسه اعتبارسنجی کارهامون خیلی بهش احتیاج دارم. پس بنابه دلایلی که بیان شد بهتره که خودم بتونم انجامش بدم. درضمن واسه حالت شبکه غیرمنظم تو سایت مث وورک کلی کد و راه حل گذاشتن واسش که اونام خیلی گیج کننده هستن و تقریبا با کار مافرق دارن چون اولا نامنظمه و راستا و لاگ و تلارانس واسش تعریف میشه و کلی دردسر داره و بعدشم داده هاش خیلی خیلی کمتر ازاین حالته. تصاویر پیوست شده مربوط به نتایج و پنجره ی ورودی گزینه ی مربوط به نمودار فوق الذکر هستش.

http://s5.picofile.com/file/8128388368/variogram.jpg
http://s5.picofile.com/file/8128388442/var_on_grid.jpg

rahnema1
چهارشنبه 11 تیر 1393, 14:07 عصر
با متلب این طور میشه

a=rand(1000000,103);
numbins=10;
maxdist2=sqrt(sum((max(a(:,1:3))-min(a(:,1:3))).^2))/2.0;
indexRandom=randperm(size(a,1),1000);
a1=a(indexRandom,:);
order = nchoosek(1:size(a1,1),2);
Xi = order(:,1);
Yi = order(:,2);
d = a1(Xi,1:3) - a1(Yi,1:3);
y = sqrt(sum (d.^2,2));
lowindex=find(y<maxdist2);
h=[y(lowindex) , ((a1(Xi(lowindex),4:end)-a1(Yi(lowindex),4:end)).^2)/2.0];
index=ceil(h(:,1)*(numbins/maxdist2));
variogram=zeros(numbins,size(a1,2)-2);
for i=1:size(variogram,2)
variogram(:,i)=accumarray(index,h(:,i),[],@mean);
end

hamidkor
چهارشنبه 11 تیر 1393, 23:55 عصر
دوست عزیز ضمن تشکر مجدد بابت زحمات بی دریغ شما باید دوباره بهتون تبریک بگم بابت شاهکار مجددتون :تشویق:

نمونه گیری تصادفی ایده ی خیلی خیلی خوبی بود و در واقع فیت این کار بودش:تشویق::تشویق:

کد آخری رو روی داده های اصلی امتحان کردم و با وجود مواجه با خطای زیر نتایجش خیلی بفکر انداختم:متفکر:

Error using accumarray
First input SUBS must contain positive integer subscripts.

درحقیقت واسه واریوگرام غیر منظم اومده میشه که گام و تعداد گام ها تعریف می شه و بعد ادامه مراحل که من نیز اینکار رو انجام دادم و همچنین نصف محدوده هم با اندازه گیری حداکثر فاصله ی ممکنه خیلی مطابقت نداشت و حداکثر رنج محدوده ی مورد بررسی را انتخاب و نصف آنرا به عنوان محدودیت اول انتخاب و محدودیت بعدی تعداد گام ها و خود گام ها بودش که کلا بصورت کد زیر از آب دراومد و نتایج هم در تصویر اتچ شده ارائه شده است. (البته نتایج نسبت به خروجی نرم افزارهای مخصوص اینکار دارای پراکندگی محسوسی است ولی خیلی خیلی امیدوار کننده است).


maxdist2=max(max(a(:,1:3))-min(a(:,1:3)))/2.0;
indexRandom=randperm(size(a,1),1000);
a1=a(indexRandom,:);
order = nchoosek(1:size(a1,1),2);
Xi = order(:,1);
Yi = order(:,2);
d = a1(Xi,1:3) - a1(Yi,1:3);
y = sqrt(sum (d.^2,2));
lowindex=find(y<maxdist2 & y==10 | y==100 | y==200 | y==300 | y==400);
h=[y(lowindex) , ((a1(Xi(lowindex),4:end)-a1(Yi(lowindex),4:end)).^2)];
u=unique(h(:,1));
c=zeros(size(u,1),size(h,2));
for i=1:numel(u)
c(i,:)=[u(i,1) sum(h(h(:,1)==u(i,1),2:end),1)];
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%
countOfu = hist(h(:,1),u);
b = [u countOfu'];
%%%%%%%%%%%%%%%%%%%%%%
g=zeros(size(c));
for i=1:size(b,1)
g(i,:)=[c(i,1) (c(i,2:end))./(2*b(i,2))];
end
plot(g(:,1),g(:,2:end))
set(gca,'XTick',g(:,1))
title('Realizations Variogram');
xlabel('Distance');
ylabel('Variogram Value');


120727

اما دو نکته:

اولا میخواستم ازخدمتتون سوال کنم که چطورمیشه که نمودار خطی رو تومتلب به صورت پلی نومیال درآورد؟ مانند تصویر سمت چپی که تو اکسل بصورت اسکاتر پلات ایجاد شده است.

دوم اینکه بنظر کارشناسانه ی شما نکته ای وجود دارد که بتوان درراستای بهبود کد بالا بکاربست یا نه؟ هرچند زمان برای 1000 انتخاب تصادفی قابل قبول بود وافزایش این مقدار هم باعث ازدیاد بیش از حد خطوط نمودار می گردد.

rahnema1
پنج شنبه 12 تیر 1393, 05:08 صبح
اگه index را به این صورت اصلاح کنید اون کدی که گذاشتم فکر کنم خطا نمیده:

index=floor(h(:,1)*(numbins/maxdist2))+1;

حتما لازم نیست درون یابی کنید. تعداد numbins را زیاد کنید مثلا 50 تا که شکلش به صورت نرم بشه
اگه numbins را زیاد کنید بهتره توی randperm هم نمونه ها را زیاد کنید مثلا بذارید 5000 تا
روش شما که y==10 | y==100 | y==200 گذاشتید فکر کنم باید یک محدوده فاصله را در بر بگیره نه یک فاصله خاص مثلا بگید بین 0 تا 100 و کلاس بعدی بین 100 تا 200 و ... که در این فاصله میانگین محاسبه بشه که من این کار را کردم

hamidkor
پنج شنبه 12 تیر 1393, 15:44 عصر
دوست عزیز تشکر مجدد، با اصلاحیه جدید جواب دادش.

درارتباط با کد خودمم باید عرض کنم خدمتتان که منظور از پلی نمیال درونیابی به معنی مستقیم نبودش، منظور تنها به قول شما نرم تر کردن نمودارها بودش.

اماقسمت هایی که اضاف شده ایده اش همانطورکه در بالا نیزاشاره شد ازاینجاست که واسه محاسبه واریوگرام معمولا گام (log) و تعداد گام درنظرگرفته می شود. و بدین خاطر از مجموعه ی فواصل محاسبه شده (h) فواصلی را که اولا کمتر از نصف حداکثر محدوده هستند (maxdist2=max(max(a(:,1:3))-min(a(:,1:3)))/2.0;) را جدا و از این بین نیز آن فواصلی که مقادیرشان 10 یا 100 یا 200 یا 300 یا 400 می باشند تعیین و واریوگرام متناظر با آنها محاسبه گردید. و در انتها ماتریس نهایی یک ماتریس با ابعاد 5 در 101 حاصل گردید.

درضمن در کد شما اگر بشه که در انتها ماتریس نهایی بجای 1 الی 10 در ستون اول فاصله ها را داشته باشد برای ایجاد نمودار خیلی مناسب خواهد بود. یعنی اگر مشخص شود که هر کلاس چقد مسافت را پوشش می دهد و کل نمودارها کلا تا چه مسافتی را پوشش داده اند خیلی مناسب خواهد بود.

یکی دو نکته دیگر که قبلا با اجازتون رو کد شما امتحان کرده بودم و زمان اجرا تاحدودی کاهش یافته بود و میخواستم نظر شما را بدانم؟ این موارد بدین صورت می باشند:


y = sum (d.^2,2);

lowindex=y<maxdist2;

h=[sqrt(y(lowindex)) , ((a1(Xi(lowindex),4:end)-a1(Yi(lowindex),4:end)).^2)];


بجای


y = sqrt(sum (d.^2,2));
lowindex=find(y<maxdist2);
h=[y(lowindex) , ((a1(Xi(lowindex),4:end)-a1(Yi(lowindex),4:end)).^2)];



باتشکر و تقدیم احترام

rahnema1
پنج شنبه 12 تیر 1393, 18:32 عصر
ستون اول که گذاشتم فاصله را میده مگه غیر از اینه؟

y = sum (d.^2,2
اینجا شما باید جذر بگیرید تا فاصله حساب بشه

lowindex=y<maxdist2
حاصل عبارت بالا یک بردار هست به طور y

اما حاصل :
lowindex=find(y<maxdist2
یک بردار به طول کوتاهتره

بنابراین اولی قاعدتا باید بیشتر طول بکشه که به عنوان اندکس قرار بگیره

hamidkor
پنج شنبه 12 تیر 1393, 19:36 عصر
درسته، و جذر گرفتن از تعداد کمی از h ها در این مرحله (h=[sqrt(y(lowindex)) , ((a1(Xi(lowindex),4:end)-a1(Yi(lowindex),4:end)).^2)];) انجام میگیره، بعدشم شاید گفتم که یه تست مختصر بود و خیلی هم تعیین کننده نیست بیخیال.

اما مسئله مهم درمورد کد شما اینه که آیا امکان گنجاندن خود فاصله ها بجای تعداد کلاس ها وجود داره یا نه؟ یعنی بجای 1 الی 10 تو ستون اول ماتریس نهایی بشه مثلا 10 الی 100 را گنجاند یا 50 الی 500.

نشد هم بیخیال کار تمومه تا اینجا.

rahnema1
پنج شنبه 12 تیر 1393, 20:50 عصر
در ستون اول فاصله قرار داره نه اعداد 1 تا 10
مثلا اگه numbins برابر 10 باشه و دامنه فاصله هم بین 0 تا 500 باشه این دامنه به ده قسمت مساوی تقسیم میشه و در این ده قسمت میانگین گرفته میشه
0-50
50-100
100-150
150-200
....
در ضمن اون جذر را هم که بعدا گرفتید عالی بود